diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index d21652c29952708d55ae33f1d0e3e18a160607a1..f1fbb16cb4fef4246f4d2e652f0e7815695d583c 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -76,8 +76,9 @@
  SCALAR managedForestAreaMaxChangeRate;
  SCALAR forestEstablishmentCost;
  SCALAR naturalAreaValue;
- SCALAR managedWoodHarvest;
- SCALAR vegClearingCost cost of clearing vegetation $1000 per tC; 
+ SCALAR woodFromRotation;
+ SCALAR vegClearingCost cost of clearing vegetation $1000 per tC;
+ SCALAR discountRate;
 
 *$gdxin "C:\Users\Bart\Documents\PhD\GAMS testing area\_gams_java_gdb1.gdx"
 $gdxin %gdxincname%
@@ -88,7 +89,7 @@ $load fertParam, irrigParam, otherIParam, exportPrices, importPrices, maxNetImpo
 $load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate
 $load previousLandCoverArea, carbonFluxRate, woodYieldLuc, carbonPrice, woodPrice
 $load conversionCost, woodYieldRota, infrastructureCost, forestRotationPeriod, woodMaxNetImport, woodMinNetImport, managedForestAreaMaxChangeRate, woodDemand
-$load forestEstablishmentCost, naturalAreaValue, managedWoodHarvest, vegClearingCost
+$load forestEstablishmentCost, naturalAreaValue, woodFromRotation, vegClearingCost, discountRate
 $gdxin
 
  SCALAR delta "use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf" / 0.00000000001 /;
@@ -98,21 +99,11 @@ $gdxin
 
  previousCropArea(crop_less_pasture, location) = previousCropArea(crop_less_pasture, location) * (1.0 - unhandledCropRate);
  
-* PARAMETER marketPrices(market_crop);
-* marketPrices('wheat') = exportPrices('wheat');
-* marketPrices('maize') = exportPrices('maize');
-* marketPrices('rice') = exportPrices('rice');
-* marketPrices('oilcrops') = exportPrices('oilcrops');
-* marketPrices('pulses') = exportPrices('pulses');
-* marketPrices('starchyRoots') = exportPrices('starchyRoots');
-* marketPrices('fruitveg') = exportPrices('fruitveg');
-* marketPrices('sugar') = exportPrices('sugar');
-* marketPrices('energycrops') = exportPrices('energycrops');
  
  PARAMETER maxExport(import_crop);
  PARAMETER minExport(import_crop);
- maxExport(import_crop) = max(-minNetImport(import_crop), 0);
- minExport(import_crop) = max(-maxNetImport(import_crop), 0);
+* maxExport(import_crop) = max(-minNetImport(import_crop), 0);
+* minExport(import_crop) = max(-maxNetImport(import_crop), 0);
 
  PARAMETER cropDM(crop)  kg DM per kg of feed the conversion from feed to meat is done in the R animal product index.  Pasture is in DM terms already
           /   wheat         0.87
@@ -143,7 +134,8 @@ $gdxin
  baseCost('pasture') = 0.04;
  otherIntCost('pasture') = 0.8 + otherICost;
  
- SCALAR discountRate / 0.08 /;
+ maxNetImport(import_crop) = min(maxNetImport(import_crop), demand(import_crop));
+ minNetImport(import_crop) = min(minNetImport(import_crop), demand(import_crop));
 
  VARIABLES
        cropArea(crop, location)           total area for crops - Mha
@@ -293,7 +285,7 @@ $gdxin
                                                           
  WOOD_HARVEST_LUC_CALC(location) .. woodHarvestLuc(location) =E= sum((forested, land_cover_after), landCoverChange(forested, land_cover_after, location) * woodYieldLuc(forested, land_cover_after, location));
 
- WOOD_SUPPLY_CALC .. woodSupply =E= managedWoodHarvest + sum(location, woodHarvestNat(location) + woodHarvestLuc(location));
+ WOOD_SUPPLY_CALC .. woodSupply =E= woodFromRotation + sum(location, woodHarvestNat(location) + woodHarvestLuc(location));
  
  WOOD_MARKET_CONSTRAINT .. woodSold =E= woodDemand + woodExported - woodImported;
  
@@ -390,7 +382,7 @@ $gdxin
  netCarbonFlux = SUM(location, carbonFlux.L(location));
  netWoodImport = woodImported.L - woodExported.L;
 * netWoodImport = 0;
- totalWoodHarvest = sum(location, woodHarvestNat.L(location) + woodHarvestLuc.L(location)) + managedWoodHarvest;
+ totalWoodHarvest = sum(location, woodHarvestNat.L(location) + woodHarvestLuc.L(location)) + woodFromRotation;
 
  Scalar totalCostsLessLU;
 
diff --git a/debug_config.properties b/debug_config.properties
index abec5b0e51e10bf4c840cbab38d160aa0da71c1b..86fd372635c363bdc216c1103f960ea646c2bb1f 100644
--- a/debug_config.properties
+++ b/debug_config.properties
@@ -7,7 +7,7 @@ OUTPUT_DIR = C:/Users/Bart/git/plumv2/output
 BASE_YEAR=2010
 START_TIMESTEP=0
 TIMESTEP_SIZE=1
-END_TIMESTEP=10
+END_TIMESTEP=5
 
 IS_CALIBRATION_RUN=true
 END_FIRST_STAGE_CALIBRATION=0
@@ -20,7 +20,7 @@ DEBUG_LIMIT_COUNTRIES=true
 DEBUG_COUNTRY_NAME=Central America
 GAMS_COUNTRY_TO_SAVE=Central America
 
-INIT_WOOD_PRICE=0.05
+INIT_WOOD_PRICE=0.4
 INIT_WOOD_STOCK=1000
 INIT_CARBON_PRICE=0.0
 FOREST_ROTATION_PERIOD=30
@@ -29,5 +29,6 @@ MANAGED_FOREST_INCREASE_COST=0.2
 MANAGED_FOREST_DECREASE_COST=0.2
 FOREST_MANAGEMENT_COST=0.2
 FOREST_MAX_CHANGE=0.02
-LAND_CHANGE_COST=0.7
-NATURAL_AREA_VALUE=0.01
\ No newline at end of file
+LAND_CHANGE_COST=0.2
+NATURAL_AREA_VALUE=0.01
+VEGETATION_CLEARING_COST=0.25
\ No newline at end of file
diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java
index e2ea12e931239ad7d7164351d52516edb8efa4c0..0a7fa9c0b3044861b8b87110a95f972b516e8754 100755
--- a/src/ac/ed/lurg/ModelConfig.java
+++ b/src/ac/ed/lurg/ModelConfig.java
@@ -451,4 +451,5 @@ public class ModelConfig {
 	public static final int FOREST_INIT_AGE = getIntProperty("FOREST_INIT_AGE", 60); 
 	public static final double DISCOUNT_RATE = getDoubleProperty("DISCOUNT_RATE", 0.08);
 	public static final double VEGETATION_CLEARING_COST = getDoubleProperty("VEGETATION_CLEARING_COST", 0.005); //$1000/tC
+	public static final boolean CONVERSION_COSTS_FROM_FILE = getBooleanProperty("CONVERSION_COSTS_FROM_FILE", false);
 }
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 65c58f4dbafd6763c1a656b22727366d6845c601..17d11761ff51d0b755e037aaab2c4009f5a27fe2 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -11,6 +11,8 @@ import java.util.Collection;
 import java.util.Map;
 import java.util.Map.Entry;
 
+import ac.ed.lurg.carbon.CarbonFluxRasterSet;
+import ac.ed.lurg.carbon.CarbonFluxReader;
 import ac.ed.lurg.country.AbstractCountryAgent;
 import ac.ed.lurg.country.AnimalRateManager;
 import ac.ed.lurg.country.CompositeCountry;
@@ -31,8 +33,6 @@ import ac.ed.lurg.forestry.ForestGrowthItem;
 import ac.ed.lurg.forestry.ForestGrowthRasterSet;
 import ac.ed.lurg.forestry.NaturalForestGrowthDataReader;
 import ac.ed.lurg.forestry.TimberForestGrowthDataReader;
-import ac.ed.lurg.landuse.CarbonFluxRasterSet;
-import ac.ed.lurg.landuse.CarbonFluxReader;
 import ac.ed.lurg.landuse.ConversionCostReader;
 import ac.ed.lurg.landuse.CropUsageData;
 import ac.ed.lurg.landuse.CropUsageReader;
@@ -42,6 +42,7 @@ import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.IrrigationMaxAmountReader;
 import ac.ed.lurg.landuse.IrrigationRasterSet;
 import ac.ed.lurg.landuse.IrrigiationCostReader;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.LandCoverItem;
 import ac.ed.lurg.landuse.LandCoverReader;
 import ac.ed.lurg.landuse.LandUseItem;
@@ -59,6 +60,7 @@ import ac.ed.lurg.types.CropType;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.types.WoodType;
 import ac.ed.lurg.utils.DoubleMap;
+import ac.ed.lurg.utils.DualKey;
 import ac.ed.lurg.utils.FileWriterHelper;
 import ac.ed.lurg.utils.LogWriter;
 import ac.ed.lurg.yield.LPJYieldResponseMapReader;
@@ -165,7 +167,8 @@ public class ModelMain {
 		
 		getForestGrowthData(timestep);
 		
-		DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts = new ConversionCostReader().read();
+		Map<LccKey, Double> conversionCosts;
+		conversionCosts = (ModelConfig.CONVERSION_COSTS_FROM_FILE) ? new ConversionCostReader().read() : new ConversionCostReader().calcFromConfig();
 		
 		countryAgents.determineProductionForAll(timestep, yieldSurfaces, currentIrrigationData, gen2Increase, currentCarbonFluxData, 
 				currentWoodYieldData, conversionCosts, carbonDemandIncrease, forestGrowthRaster);
diff --git a/src/ac/ed/lurg/carbon/AgeFunction.java b/src/ac/ed/lurg/carbon/AgeFunction.java
new file mode 100644
index 0000000000000000000000000000000000000000..dd9320356a69b3d6340206d1a23d6c510c00ed34
--- /dev/null
+++ b/src/ac/ed/lurg/carbon/AgeFunction.java
@@ -0,0 +1,40 @@
+package ac.ed.lurg.carbon;
+
+import java.io.Serializable;
+import java.util.HashMap;
+import java.util.Map;
+
+public class AgeFunction implements Serializable {
+	private static final long serialVersionUID = -6068915999792275002L;
+	
+	private Map<Integer, Double> data = new HashMap<Integer, Double>();
+	private int domainLimit = 0;
+	
+	public void setData(int age, double value) {
+		data.put(age, value);
+		domainLimit = Math.max(domainLimit, age);
+	}
+	
+	public double getValue(int age) {
+		if (data.containsKey(age)) {
+			return data.get(age);
+		}
+		if (age >= domainLimit) {
+			return data.get(domainLimit);
+		}
+		return interpolateValue(age);
+	}
+	
+	private double interpolateValue(int age) {
+		int lowerAge = data.keySet().stream().filter(v -> v.intValue() < age).reduce((result, current) -> age - current < age - result ? current : result).get();
+		int upperAge = data.keySet().stream().filter(v -> v.intValue() > age).reduce((result, current) -> current - age < result - age ? current : result).get();
+		double gradient = (data.get(upperAge) - data.get(lowerAge)) / (upperAge - lowerAge);
+		double d = data.get(lowerAge) + (age - lowerAge) * gradient;
+		return d;
+	}
+
+	public int getDomainLimit() {
+		return domainLimit;
+	}
+	
+}
diff --git a/src/ac/ed/lurg/carbon/CarbonFluxItem.java b/src/ac/ed/lurg/carbon/CarbonFluxItem.java
new file mode 100644
index 0000000000000000000000000000000000000000..2725db46c560d68fcd9fd3b1b536528465ef72b2
--- /dev/null
+++ b/src/ac/ed/lurg/carbon/CarbonFluxItem.java
@@ -0,0 +1,77 @@
+package ac.ed.lurg.carbon;
+
+import ac.sac.raster.RasterItem;
+
+import java.io.Serializable;
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.landuse.LccKey;
+import ac.ed.lurg.types.LandCoverType;
+
+public class CarbonFluxItem implements RasterItem, Serializable {
+	private static final long serialVersionUID = 440720456140537815L;
+	
+	Map<LandCoverType, Map<LandCoverType, Double>> carbonFluxes = new HashMap<LandCoverType, Map<LandCoverType, Double>>(); 
+	Map<LccKey, AgeFunction> conversionCarbonFlux = new HashMap<LccKey, AgeFunction>();
+	Map<LandCoverType, AgeFunction> ecosystemCarbonFlux = new HashMap<LandCoverType, AgeFunction>();
+	
+	public void setConversionCarbonFlux(LandCoverType fromLc, LandCoverType toLc, int age, double cFlux) {
+		LccKey key = new LccKey(fromLc, toLc);
+		AgeFunction ageFunc = conversionCarbonFlux.get(key);
+		if (ageFunc == null) {
+			ageFunc = new AgeFunction();
+			conversionCarbonFlux.put(key, ageFunc);
+		} 
+		ageFunc.setData(age, cFlux); //TODO does this work?	
+	}
+	
+	public void setEcosystemCarbonFlux(LandCoverType lcType, int age, double cFlux) {
+		AgeFunction ageFunc = ecosystemCarbonFlux.get(lcType);
+		if (ageFunc == null) {
+			ageFunc = new AgeFunction();
+			ecosystemCarbonFlux.put(lcType, ageFunc);
+		} 
+		ageFunc.setData(age, cFlux);
+	}
+	
+	public double getConversionCarbonFlux(LccKey key, int age) {
+		return conversionCarbonFlux.get(key).getValue(age);
+	}
+	
+	public double getEcosystemCarbonFlux(LandCoverType lcType, int age) {
+		return ecosystemCarbonFlux.get(lcType).getValue(age);
+	}
+	
+	public void setCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover, double carbonFlux) {
+		
+		if (carbonFluxes.containsKey(previousLandCover)) {
+			carbonFluxes.get(previousLandCover).put(newLandCover, carbonFlux);
+		} else {
+			Map<LandCoverType, Double> cfMap = new HashMap<LandCoverType, Double>();
+			cfMap.put(newLandCover, carbonFlux);
+			carbonFluxes.put(previousLandCover, cfMap);
+		}
+	}
+	
+	public double getCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover) {
+		return carbonFluxes.get(previousLandCover).get(newLandCover);
+	}
+	
+	public boolean checkForKeys(LandCoverType previousLandCover, LandCoverType newLandCover) {
+		if (carbonFluxes.containsKey(previousLandCover)) {
+			if (carbonFluxes.get(previousLandCover).containsKey(newLandCover)) {
+				return true;
+			} else {
+				return false;
+			}
+		} else {
+			return false;
+		}
+	}
+	
+	public Map<LandCoverType, Map<LandCoverType, Double>> getCarbonFluxes() {
+		return carbonFluxes;
+	}
+
+}
diff --git a/src/ac/ed/lurg/landuse/CarbonFluxRasterSet.java b/src/ac/ed/lurg/carbon/CarbonFluxRasterSet.java
similarity index 96%
rename from src/ac/ed/lurg/landuse/CarbonFluxRasterSet.java
rename to src/ac/ed/lurg/carbon/CarbonFluxRasterSet.java
index 986bea65a2feab9b8482f100fc7c44ef9057cf20..3289872a0277b89877eaf398886432c08e5ac976 100644
--- a/src/ac/ed/lurg/landuse/CarbonFluxRasterSet.java
+++ b/src/ac/ed/lurg/carbon/CarbonFluxRasterSet.java
@@ -1,4 +1,4 @@
-package ac.ed.lurg.landuse;
+package ac.ed.lurg.carbon;
 
 import java.util.Collection;
 
diff --git a/src/ac/ed/lurg/landuse/CarbonFluxReader.java b/src/ac/ed/lurg/carbon/CarbonFluxReader.java
similarity index 99%
rename from src/ac/ed/lurg/landuse/CarbonFluxReader.java
rename to src/ac/ed/lurg/carbon/CarbonFluxReader.java
index 197ee781895709d3f5270a1cfedf22fbb8043a10..9b71255736d51437b686c907d4bfc765f3c80c64 100644
--- a/src/ac/ed/lurg/landuse/CarbonFluxReader.java
+++ b/src/ac/ed/lurg/carbon/CarbonFluxReader.java
@@ -1,4 +1,4 @@
-package ac.ed.lurg.landuse;
+package ac.ed.lurg.carbon;
 
 import java.util.Map;
 
diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java
index 3b59117dc8ef52245862974cf97370be05187a40..9266a5f18b28da23db691348050c42744d3dfeb6 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -5,22 +5,24 @@ import java.io.IOException;
 import java.nio.file.FileSystems;
 import java.nio.file.Files;
 import java.nio.file.StandardCopyOption;
+import java.util.Collection;
 import java.util.HashMap;
 import java.util.HashSet;
 import java.util.List;
 import java.util.Map;
 import java.util.Map.Entry;
 import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.carbon.CarbonFluxItem;
 import ac.ed.lurg.country.gams.GamsCountryInput;
 import ac.ed.lurg.country.gams.GamsRasterInput;
 import ac.ed.lurg.country.gams.GamsRasterOptimiser;
 import ac.ed.lurg.country.gams.GamsRasterOutput;
 import ac.ed.lurg.demand.AbstractDemandManager;
 import ac.ed.lurg.forestry.ForestGrowthItem;
-import ac.ed.lurg.landuse.CarbonFluxItem;
 import ac.ed.lurg.landuse.CropUsageData;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.WoodUsageData;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.types.CommodityType;
@@ -108,7 +110,7 @@ public class CountryAgent extends AbstractCountryAgent {
 	
 	public GamsRasterOutput determineProduction(YieldRaster countryYieldSurfaces, RasterSet<IrrigationItem> irrigData, 
 			Map<CropType, GlobalPrice> worldPrices, double globalGen2EcIncrease, RasterSet<CarbonFluxItem> carbonFluxData, 
-			RasterSet<WoodYieldItem> woodYieldData, DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts,
+			RasterSet<WoodYieldItem> woodYieldData, Map<LccKey, Double> conversionCosts,
 			double carbonDemandIncrease, GlobalPrice carbonPrice, GlobalPrice woodPrice, RasterSet<ForestGrowthItem> forestGrowthData) {
 
 		// project demand
@@ -176,7 +178,7 @@ public class CountryAgent extends AbstractCountryAgent {
 	}
 
 	private GamsRasterInput getGamsRasterInput(RasterSet<IrrigationItem> irrigData, YieldRaster countryYieldSurfaces, double gen2EcIncrease,
-			RasterSet<CarbonFluxItem> carbonFluxData, RasterSet<WoodYieldItem> woodYieldData, DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts,
+			RasterSet<CarbonFluxItem> carbonFluxData, RasterSet<WoodYieldItem> woodYieldData, Map<LccKey, Double> conversionCosts,
 			double carbonDemandIncrease) {
 		double allowedImportChange;
 
@@ -186,6 +188,17 @@ public class CountryAgent extends AbstractCountryAgent {
 		else { // normal (not the initial) time-step
 			allowedImportChange = ModelConfig.MAX_IMPORT_CHANGE;  // when running is calibration model calibrate (ModelConfig.IS_CALIBRATION_RUN==true) MAX_IMPORT_CHANGE is already 0.
 		}
+		
+		Map<CropType, Double> cropDemand = new HashMap<CropType, Double>();
+		for (Map.Entry<CommodityType, Double> entry : currentProjectedDemand.entrySet()) {
+			CommodityType commodity = entry.getKey();
+			double demand = entry.getValue();
+			Collection<CropType> crops = commodity.getCropTypes();
+			for (CropType c : crops) {
+				double minFract = (currentMinDemandFract.containsKey(commodity)) ? currentMinDemandFract.get(commodity).get(c) : 1;
+				cropDemand.put(c, demand * minFract);
+			}
+		}
 
 		Map<CropType, TradeConstraint> importConstraints = new HashMap<CropType, TradeConstraint>();
 
@@ -193,6 +206,10 @@ public class CountryAgent extends AbstractCountryAgent {
 			CropUsageData cropUsage = entry.getValue();
 			CropType crop = entry.getKey();
 			double baseTrade = cropUsage.getNetImportsExpected();
+			// Make sure max import <= demand so GAMS optimisation doesn't fail
+			if (ModelConfig.IS_CALIBRATION_RUN && currentTimestep.getYear() <= ModelConfig.END_FIRST_STAGE_CALIBRATION) {
+				baseTrade = Math.min(baseTrade, cropDemand.get(crop));
+			} 			
 			double changeUp = 0.0;
 			double changeDown = 0.0;
 
@@ -368,7 +385,7 @@ public class CountryAgent extends AbstractCountryAgent {
 		for (LandUseItem luItem : previousLandUses.values()) {
 			if (luItem == null)
 				continue;
-			luItem.growForests();
+			luItem.updateTilesAge();
 		}
 	}
 	
diff --git a/src/ac/ed/lurg/country/CountryAgentManager.java b/src/ac/ed/lurg/country/CountryAgentManager.java
index 5285d49c0232ef564c7ac04a2d3f8f817395f85e..b2be34efdd98ec7beff6c8138befe7400d988721 100644
--- a/src/ac/ed/lurg/country/CountryAgentManager.java
+++ b/src/ac/ed/lurg/country/CountryAgentManager.java
@@ -12,17 +12,18 @@ import java.util.Map;
 import ac.ed.lurg.InternationalMarket;
 import ac.ed.lurg.ModelConfig;
 import ac.ed.lurg.Timestep;
+import ac.ed.lurg.carbon.CarbonFluxItem;
+import ac.ed.lurg.carbon.CarbonFluxRasterSet;
 import ac.ed.lurg.country.crafty.CraftyCountryAgent;
 import ac.ed.lurg.country.crafty.CraftyProdManager;
 import ac.ed.lurg.demand.AbstractDemandManager;
 import ac.ed.lurg.forestry.ForestGrowthItem;
 import ac.ed.lurg.forestry.ForestGrowthRasterSet;
-import ac.ed.lurg.landuse.CarbonFluxItem;
-import ac.ed.lurg.landuse.CarbonFluxRasterSet;
 import ac.ed.lurg.landuse.CropUsageData;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.IrrigationRasterSet;
 import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.WoodUsageData;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.landuse.WoodYieldRasterSet;
@@ -109,7 +110,7 @@ public class CountryAgentManager {
 	}
 
 	public void determineProductionForAll(Timestep timestep, YieldRaster yieldSurfaces, IrrigationRasterSet currentIrrigationData, double gen2Increase,
-			CarbonFluxRasterSet currentCarbonFluxData, WoodYieldRasterSet currentWoodYieldData, DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts,
+			CarbonFluxRasterSet currentCarbonFluxData, WoodYieldRasterSet currentWoodYieldData, Map<LccKey, Double> conversionCosts,
 			double carbonDemandIncrease, ForestGrowthRasterSet forestGrowthRaster) {
 		
 		for (AbstractCountryAgent aca : countryAgents) {		
diff --git a/src/ac/ed/lurg/country/LandCoverChangeItem.java b/src/ac/ed/lurg/country/LandCoverChangeItem.java
index 1c01df59632bd44001283ded0a9c8e408990125b..5e06a36fc9b70419f5f9ecdc52a17f4400a08559 100644
--- a/src/ac/ed/lurg/country/LandCoverChangeItem.java
+++ b/src/ac/ed/lurg/country/LandCoverChangeItem.java
@@ -1,39 +1,35 @@
 package ac.ed.lurg.country;
 
+import ac.ed.lurg.types.LandCoverType;
+
 public class LandCoverChangeItem {
 	private int location;
-	private String fromLandCover;
-	private String toLandCover;
+	private LandCoverType fromLandCover;
+	private LandCoverType toLandCover;
 	private double area;
-	private double timberYield;
 	
 	public LandCoverChangeItem() {}
 	
-	public LandCoverChangeItem(int location, String fromLandCover, String toLandCover, double area, double timberYield) {
+	public LandCoverChangeItem(int location, LandCoverType fromLandCover, LandCoverType toLandCover, double area) {
 		this.location = location;
 		this.fromLandCover = fromLandCover;
 		this.toLandCover = toLandCover;
 		this.area = area;
-		this.timberYield = timberYield;
 	}
 
 	public int getLocation() {
 		return location;
 	}
 
-	public String getFromLandCover() {
+	public LandCoverType getFromLandCover() {
 		return fromLandCover;
 	}
 
-	public String getToLandCover() {
+	public LandCoverType getToLandCover() {
 		return toLandCover;
 	}
 
 	public double getArea() {
 		return area;
 	}
-
-	public double getTimberYield() {
-		return timberYield;
-	}
 }
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationInput.java b/src/ac/ed/lurg/country/gams/GamsLocationInput.java
index b89f020062dfe57c5e816634738155332780fbed..5ed7d25509f2378f49dd2bf5c7bc83e40d7dafb4 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationInput.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationInput.java
@@ -3,9 +3,10 @@ package ac.ed.lurg.country.gams;
 import java.util.Map;
 
 import ac.ed.lurg.Timestep;
-import ac.ed.lurg.landuse.CarbonFluxItem;
+import ac.ed.lurg.carbon.CarbonFluxItem;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.utils.DoubleMap;
@@ -19,14 +20,14 @@ public class GamsLocationInput {
 	private Map<Integer, ? extends IrrigationItem> irrigationCosts;
 	private Map<Integer, ? extends CarbonFluxItem> carbonFluxes;
 	private Map<Integer, ? extends WoodYieldItem> woodYields;
-	private DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts;
+	private Map<LccKey, Double> conversionCosts;
 	private GamsCountryInput countryInput;
 	private Map<Integer, Double> rotationPeriods;
 	
 	public GamsLocationInput(Timestep timestep, Map<Integer, ? extends YieldResponsesItem> yields, Map<Integer, ? extends LandUseItem> previousLandUse,
 			Map<Integer, ? extends IrrigationItem> irrigationCosts, Map<Integer, ? extends CarbonFluxItem> carbonFluxes, 
 			Map<Integer, ? extends WoodYieldItem> woodYields,
-			DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts, GamsCountryInput countryInput, Map<Integer, Double> rotationPeriods) {
+			Map<LccKey, Double> conversionCosts, GamsCountryInput countryInput, Map<Integer, Double> rotationPeriods) {
 		super();
 		this.timestep = timestep;
 		this.yields = yields;
@@ -59,7 +60,7 @@ public class GamsLocationInput {
 		return woodYields;
 	}
 	
-	public DoubleMap<LandCoverType, LandCoverType, Double> getConversionCosts() {
+	public Map<LccKey, Double> getConversionCosts() {
 		return conversionCosts;
 	}
 	
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index 2b2cac5cf996eae209cb18e6694da43511d75387..6f6c97c6c3020d2f0155759d76a6e98d6a70d192 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -22,14 +22,15 @@ import com.gams.api.GAMSWorkspace;
 import com.gams.api.GAMSWorkspaceInfo;
 
 import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.carbon.CarbonFluxItem;
 import ac.ed.lurg.country.CountryPrice;
 import ac.ed.lurg.country.LandCoverChangeItem;
 import ac.ed.lurg.country.TradeConstraint;
-import ac.ed.lurg.landuse.CarbonFluxItem;
 import ac.ed.lurg.landuse.CropUsageData;
 import ac.ed.lurg.landuse.Intensity;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.WoodUsageData;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.types.CommodityType;
@@ -409,50 +410,17 @@ public class GamsLocationOptimiser {
 		
 		// Land cover conversion cost
 		GAMSParameter conversionCostP = inDB.addParameter("conversionCost", 2);
-		/*
-		DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts = inputData.getConversionCosts();
-		
-		*/
-		
-		DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts = new DoubleMap<LandCoverType, LandCoverType, Double>();
-		
-		conversionCosts.put(LandCoverType.CROPLAND, LandCoverType.PASTURE, ModelConfig.CROP_DECREASE_COST + ModelConfig.PASTURE_INCREASE_COST);
-		conversionCosts.put(LandCoverType.CROPLAND, LandCoverType.NATURAL, ModelConfig.CROP_DECREASE_COST);
-		conversionCosts.put(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST, ModelConfig.CROP_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
-		conversionCosts.put(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST, ModelConfig.CROP_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
-		
-		conversionCosts.put(LandCoverType.PASTURE, LandCoverType.CROPLAND, ModelConfig.PASTURE_DECREASE_COST + ModelConfig.CROP_INCREASE_COST);
-		conversionCosts.put(LandCoverType.PASTURE, LandCoverType.NATURAL, ModelConfig.PASTURE_DECREASE_COST);
-		conversionCosts.put(LandCoverType.PASTURE, LandCoverType.TIMBER_FOREST, ModelConfig.PASTURE_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
-		conversionCosts.put(LandCoverType.PASTURE, LandCoverType.CARBON_FOREST, ModelConfig.PASTURE_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
-		
-		conversionCosts.put(LandCoverType.NATURAL, LandCoverType.CROPLAND, ModelConfig.AGRI_EXPANSION_COST_BASE + ModelConfig.CROP_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
-		conversionCosts.put(LandCoverType.NATURAL, LandCoverType.PASTURE, ModelConfig.AGRI_EXPANSION_COST_BASE + ModelConfig.PASTURE_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
-		conversionCosts.put(LandCoverType.NATURAL, LandCoverType.TIMBER_FOREST, ModelConfig.MANAGED_FOREST_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
-		conversionCosts.put(LandCoverType.NATURAL, LandCoverType.CARBON_FOREST, ModelConfig.MANAGED_FOREST_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
-		
-		conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.CROPLAND, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
-		conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.PASTURE, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.PASTURE_INCREASE_COST);
-		conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.NATURAL, ModelConfig.MANAGED_FOREST_DECREASE_COST);
-		conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST, 0.0);
-		conversionCosts.put(LandCoverType.TIMBER_FOREST, LandCoverType.CARBON_FOREST, ModelConfig.FOREST_MANAGEMENT_COST);
-		
-		conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.CROPLAND, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
-		conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.PASTURE, ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.PASTURE_INCREASE_COST);
-		conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.NATURAL, ModelConfig.MANAGED_FOREST_DECREASE_COST);
-		conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.TIMBER_FOREST, ModelConfig.FOREST_MANAGEMENT_COST);
-		conversionCosts.put(LandCoverType.CARBON_FOREST, LandCoverType.CARBON_FOREST, ModelConfig.FOREST_MANAGEMENT_COST);
-		
-		for (Entry<LandCoverType, Map<LandCoverType, Double>> fromMap : conversionCosts.getMap().entrySet()) {
-			String fromName = fromMap.getKey().getName();
-			for (Entry<LandCoverType, Double> toMap : fromMap.getValue().entrySet()) {
-				String toName = toMap.getKey().getName();
-				double cost = toMap.getValue();
-				Vector<String> v = new Vector<String>();
-				v.add(fromName);
-				v.add(toName);
-				setGamsParamValue(conversionCostP.addRecord(v), cost, 5);
-			}
+		Map<LccKey, Double> conversionCosts = inputData.getConversionCosts();
+		
+		for (Map.Entry<LccKey, Double> entry : conversionCosts.entrySet()) {
+			LccKey key = entry.getKey();
+			String fromName = key.getFromLc().getName();
+			String toName = key.getToLc().getName();
+			double cost = entry.getValue();
+			Vector<String> v = new Vector<String>();
+			v.add(fromName);
+			v.add(toName);
+			setGamsParamValue(conversionCostP.addRecord(v), cost, 5);
 		}
 		
 		GAMSParameter forestRotationP = inDB.addParameter("forestRotationPeriod", 1);
@@ -468,8 +436,9 @@ public class GamsLocationOptimiser {
 		addScalar(inDB, "forestEstablishmentCost", ModelConfig.FOREST_ESTABLISHMENT_COST, 5);
 		addScalar(inDB, "naturalAreaValue", ModelConfig.NATURAL_AREA_VALUE, 5);
 		double managedWoodHarvest = countryInput.getWoodHarvestFromRotation();
-		addScalar(inDB, "managedWoodHarvest", managedWoodHarvest, 5);
+		addScalar(inDB, "woodFromRotation", managedWoodHarvest, 5);
 		addScalar(inDB, "vegClearingCost", ModelConfig.VEGETATION_CLEARING_COST, 5);
+		addScalar(inDB, "discountRate", ModelConfig.DISCOUNT_RATE, 5);
 	}
 
 	private void addScalar(GAMSDatabase gamsDb, String recordName, double val, int places) {
@@ -581,7 +550,7 @@ public class GamsLocationOptimiser {
 				totalCropArea += area;
 			}
 			
-			landUseItem.setCropFraction(cropType, croplandArea > 0 ? area/croplandArea : 0);
+			landUseItem.setCropFraction(cropType, croplandArea > 0 ? area/croplandArea : 0); // TODO gives wrong pasture numbers
 			
 		}
 
@@ -600,34 +569,24 @@ public class GamsLocationOptimiser {
 			
 		// Land cover change
 		GAMSVariable varLandCoverChange = outDB.getVariable("landCoverChange");
-		
-		TripleMap<Integer, LandCoverType, LandCoverType, Double> landCoverChanges = new TripleMap<Integer, LandCoverType, LandCoverType, Double>();
-		ArrayList<LandCoverChangeItem> gamsLandCoverChanges = new ArrayList<LandCoverChangeItem>();
+		Map<Integer, ArrayList<LandCoverChangeItem>> landCoverChanges = new HashMap<Integer, ArrayList<LandCoverChangeItem>>();
 		
 		for (GAMSVariableRecord rec : varLandCoverChange) {
-			String fromGamsLc = rec.getKeys()[0];
-			String toGamsLc = rec.getKeys()[1];
+			String fromLcStr = rec.getKeys()[0];
+			String toLcStr = rec.getKeys()[1];
 			String locationName = rec.getKeys()[2];
 			int locId = Integer.parseInt(locationName);
 			double change = rec.getLevel();
 			
 			// Skip if no change
-			if (change == 0) 
+			if (change == 0 || fromLcStr.equals(toLcStr)) 
 				continue;
 
-			// Data for updating land cover
-			String fromLcStr = (fromGamsLc.equals("naturalNew")) ? "natural" : fromGamsLc;
-			String toLcStr = (toGamsLc.equals("naturalNew")) ? "natural" : toGamsLc;
 			LandCoverType fromLc = LandCoverType.getForName(fromLcStr);
-			LandCoverType toLc = LandCoverType.getForName(toLcStr);		
-			
-			if (!fromLc.equals(toLc)) { //Important as don't want to include unchanged land cover
-				landCoverChanges.put(locId, fromLc, toLc, change);
-			}
+			LandCoverType toLc = LandCoverType.getForName(toLcStr);	
 			
-			// Data for updating natural areas and forest
-			double potentialTimberYield = (toLc.equals(LandCoverType.TIMBER_FOREST)) ? inputData.getWoodYields().get(locId).getWoodYieldRota(fromLc, toLc) : 0;
-			gamsLandCoverChanges.add(new LandCoverChangeItem(locId, fromGamsLc, toGamsLc, change, potentialTimberYield));
+			ArrayList<LandCoverChangeItem> changesList = landCoverChanges.computeIfAbsent(locId, k -> new ArrayList<LandCoverChangeItem>());
+			changesList.add(new LandCoverChangeItem(locId, fromLc, toLc, change));
 		}
 
 		// Carbon flux
@@ -661,7 +620,7 @@ public class GamsLocationOptimiser {
 			naturalForestCleared.put(locId, areaCleared);
 		}
 
-		GamsLocationOutput results = new GamsLocationOutput(modelStatus, landUses, cropUsageData, landCoverChanges, gamsLandCoverChanges, netCarbonFlux, newWoodUsageMap, naturalForestCleared);
+		GamsLocationOutput results = new GamsLocationOutput(modelStatus, landUses, cropUsageData, landCoverChanges, netCarbonFlux, newWoodUsageMap, naturalForestCleared);
 		return results;
 	}
 	
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOutput.java b/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
index 1d1b177f4f987e77537538095efd1cc4aea0be17..edce072f7f69549863c58b4e99b2145e888b2eaf 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOutput.java
@@ -19,8 +19,7 @@ public class GamsLocationOutput {
 	
 	Map<Integer, LandUseItem> landUses;  // data mapped from id (not raster)
 	private Map<CropType, CropUsageData> cropUsageData;
-	private TripleMap<Integer, LandCoverType, LandCoverType, Double> landCoverChanges;
-	ArrayList<LandCoverChangeItem> gamsLandCoverChanges;
+	Map<Integer, ArrayList<LandCoverChangeItem>> landCoverChanges;
 	private double netCarbonFlux;
 	private Map<WoodType, WoodUsageData> woodUsageData;
 	private Map<Integer, Double> naturalForestCleared;
@@ -28,16 +27,14 @@ public class GamsLocationOutput {
 	public GamsLocationOutput(ModelStat status, 
 			Map<Integer, LandUseItem> landUses, 
 			Map<CropType, CropUsageData> cropUsageData,
-			TripleMap<Integer, LandCoverType, LandCoverType, Double> landCoverChange,
-			ArrayList<LandCoverChangeItem> gamsLandCoverChanges,
+			Map<Integer, ArrayList<LandCoverChangeItem>> landCoverChanges,
 			double netCarbonFlux, Map<WoodType, WoodUsageData> woodUsageData,
 			Map<Integer, Double> naturalForestCleared) {
 		super();
 		this.status = status;
 		this.landUses = landUses;
 		this.cropUsageData = cropUsageData;
-		this.landCoverChanges = landCoverChange;
-		this.gamsLandCoverChanges = gamsLandCoverChanges;
+		this.landCoverChanges = landCoverChanges;
 		this.netCarbonFlux = netCarbonFlux;
 		this.woodUsageData = woodUsageData;
 		this.naturalForestCleared = naturalForestCleared;
@@ -54,14 +51,10 @@ public class GamsLocationOutput {
 		return cropUsageData;
 	}
 	
-	public TripleMap<Integer, LandCoverType, LandCoverType, Double> getLandCoverChanges() {
+	public Map<Integer, ArrayList<LandCoverChangeItem>> getLandCoverChanges() {
 		return landCoverChanges;
 	}
 
-	public ArrayList<LandCoverChangeItem> getGamsLandCoverChanges() {
-		return gamsLandCoverChanges;
-	}
-
 	public double getNetCarbonFlux() {
 		return netCarbonFlux;
 	}
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterInput.java b/src/ac/ed/lurg/country/gams/GamsRasterInput.java
index 3bf7eb73649cae1ece1e26751453729aa6f20fd7..07b9d1de3804fea288271ccd9936649593358fae 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterInput.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterInput.java
@@ -1,11 +1,14 @@
 package ac.ed.lurg.country.gams;
 
+import java.util.Map;
+
 import ac.ed.lurg.Timestep;
+import ac.ed.lurg.carbon.CarbonFluxItem;
 import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.landuse.LccKey;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.utils.DoubleMap;
-import ac.ed.lurg.landuse.CarbonFluxItem;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.yield.YieldRaster;
 import ac.sac.raster.RasterSet;
@@ -18,12 +21,12 @@ public class GamsRasterInput {
 	private RasterSet<IrrigationItem> irrigationCost;
 	private RasterSet<CarbonFluxItem> carbonFluxes;
 	private RasterSet<WoodYieldItem> woodYields;
-	private DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts;
+	private Map<LccKey, Double> conversionCosts;
 	private GamsCountryInput countryInput;
 
 	public GamsRasterInput(Timestep timestep, YieldRaster yields, RasterSet<LandUseItem> previousLandsUses, RasterSet<IrrigationItem> irrigationCost, 
 			RasterSet<CarbonFluxItem> carbonFluxes, RasterSet<WoodYieldItem> woodYields,
-			DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts, GamsCountryInput countryInput) {
+			Map<LccKey, Double> conversionCosts, GamsCountryInput countryInput) {
 		super();
 		this.timestep = timestep;
 		this.yields = yields;
@@ -55,7 +58,7 @@ public class GamsRasterInput {
 		return woodYields;
 	}
 	
-	public DoubleMap<LandCoverType, LandCoverType, Double> getConversionCosts() {
+	public Map<LccKey, Double> getConversionCosts() {
 		return conversionCosts;
 	}
 
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
index 43dda5d371c49a4be79f301b2de57c533a47ff29..e19783d9f5a7739fc6a583ed5b2a1bdd89d43d60 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
@@ -1,12 +1,14 @@
 package ac.ed.lurg.country.gams;
 
+import java.util.ArrayList;
 import java.util.HashMap;
 import java.util.HashSet;
 import java.util.Map;
 import java.util.Map.Entry;
 import java.util.Set;
 
-import ac.ed.lurg.landuse.CarbonFluxItem;
+import ac.ed.lurg.carbon.CarbonFluxItem;
+import ac.ed.lurg.country.LandCoverChangeItem;
 import ac.ed.lurg.landuse.Intensity;
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
@@ -58,7 +60,7 @@ public class GamsRasterOptimiser {
 			baseTimberYield.put(locId, gamsInput.getWoodYields().get(locId).getWoodYieldRota(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST));
 		}
 
-		return new GamsRasterOutput(gamsOutput.getStatus(), newIntensityRaster, gamsOutput.getCommoditiesData(), gamsOutput.getGamsLandCoverChanges(), 
+		return new GamsRasterOutput(gamsOutput.getStatus(), newIntensityRaster, gamsOutput.getCommoditiesData(),
 				baseTimberYield, gamsOutput.getNetCarbonFlux(), gamsOutput.getWoodUsageData());
 	}
 
@@ -118,17 +120,16 @@ public class GamsRasterOptimiser {
 			
 			RasterSet<LandUseItem> landUseItemsForLocation = newLandUseRaster.createSubsetForKeys(keys);
 			
-			DoubleMap<LandCoverType, LandCoverType, Double> landCoverChange = gamsOutput.getLandCoverChanges().getInnerMap(locId);
+			ArrayList<LandCoverChangeItem> landCoverChanges = gamsOutput.getLandCoverChanges().get(locId);
 			
-			if (landCoverChange != null) {	
+			if (landCoverChanges != null) {	
 				// Allocate land cover change to grid cells
-				for (Entry<LandCoverType, Map<LandCoverType, Double>> fromMap: landCoverChange.getMap().entrySet()) {
-					LandCoverType fromLC = fromMap.getKey();
-					for (Entry<LandCoverType, Double> toMap: fromMap.getValue().entrySet()) {
-						LandCoverType toLC = toMap.getKey();
-						double change = toMap.getValue();
-						allocAllLandCropsTop(newLandUseRaster, keys, toLC, fromLC, change, locId);
-					}
+				for (LandCoverChangeItem lccItem : landCoverChanges) {
+					LandCoverType fromLC = lccItem.getFromLandCover();
+					LandCoverType toLC = lccItem.getToLandCover();
+					double change = lccItem.getArea();
+					allocAllLandCropsTop(newLandUseRaster, keys, toLC, fromLC, change, locId);
+					
 				}
 			}
 			
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOutput.java b/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
index 0add17f98234b42295ce9211b31726c595ef873c..2204badb514fb1e20c5e4c462665ca219d39c1dd 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOutput.java
@@ -18,7 +18,6 @@ public class GamsRasterOutput {
 	private ModelStat status;
 	private RasterSet<LandUseItem> landUses;
 	private Map<CropType, CropUsageData> cropUsageData;
-	private ArrayList<LandCoverChangeItem> gamsLandCoverChanges;
 	private Map<Integer, Double> baseTimberYield;
 	private double netCarbonFlux;
 	private Map<WoodType, WoodUsageData> woodUsageData;
@@ -37,10 +36,9 @@ public class GamsRasterOutput {
 	}
 */
 	public GamsRasterOutput(ModelStat status, RasterSet<LandUseItem> intensityRaster, Map<CropType, CropUsageData> cropUsageData,
-			ArrayList<LandCoverChangeItem> gamsLandCoverChanges, Map<Integer, Double> baseTimberYield, double netCarbonFlux, Map<WoodType, WoodUsageData> woodUsageData) {
+			Map<Integer, Double> baseTimberYield, double netCarbonFlux, Map<WoodType, WoodUsageData> woodUsageData) {
 		this(intensityRaster, cropUsageData, woodUsageData);
 		this.status = status;
-		this.gamsLandCoverChanges = gamsLandCoverChanges;
 		this.baseTimberYield = baseTimberYield;
 		this.netCarbonFlux = netCarbonFlux;
 	}
@@ -56,10 +54,6 @@ public class GamsRasterOutput {
 	public Map<CropType, CropUsageData> getCropUsageData() {
 		return cropUsageData;
 	}
-	
-	public ArrayList<LandCoverChangeItem> getGamsLandCoverChanges() {
-		return gamsLandCoverChanges;
-	}
 
 	public Map<Integer, Double> getBaseTimberYield() {
 		return baseTimberYield;
diff --git a/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java
index 8b65ec9df73137c13146c4edaa04912c5653547a..9e55bbc5f85d046ba5e1aecf18c17ccd4325a990 100644
--- a/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java
+++ b/src/ac/ed/lurg/forestry/CarbonForestGrowthDataReader.java
@@ -16,27 +16,27 @@ public class CarbonForestGrowthDataReader extends AbstractTabularRasterReader<Fo
 	}
 
 	protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 0, 0.0);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 5,  getValueForCol(rowValues, "5") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.CARBON_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 0, 0.0);
+		item.setYield(LandCoverType.CARBON_FOREST, 5,  getValueForCol(rowValues, "5") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
+		item.setYield(LandCoverType.CARBON_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
 
                                                                         
 	}
diff --git a/src/ac/ed/lurg/forestry/ForestGrowthItem.java b/src/ac/ed/lurg/forestry/ForestGrowthItem.java
index 0fc926a6e12ac2e5b7d91815ccac63c5762e8423..6a8280585019c12893990b3ced21f52dea4f1035 100644
--- a/src/ac/ed/lurg/forestry/ForestGrowthItem.java
+++ b/src/ac/ed/lurg/forestry/ForestGrowthItem.java
@@ -2,8 +2,10 @@ package ac.ed.lurg.forestry;
 
 import java.io.Serializable;
 import java.util.Collections;
+import java.util.HashMap;
 import java.util.Map;
 
+import ac.ed.lurg.carbon.AgeFunction;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.utils.DoubleMap;
 import ac.ed.lurg.utils.Interpolator;
@@ -12,43 +14,34 @@ import ac.sac.raster.RasterItem;
 public class ForestGrowthItem implements RasterItem, Serializable {
 	private static final long serialVersionUID = -4822328490070047362L;
 	
-	private DoubleMap<LandCoverType, Integer, Double> forestGrowthFunction = new DoubleMap<LandCoverType, Integer, Double>();
+	private Map<LandCoverType, AgeFunction> woodYields = new HashMap<LandCoverType, AgeFunction>();
 	
-	public void setForestGrowthFunction(LandCoverType forestType, int time, double deltaC) {
-		forestGrowthFunction.put(forestType, time, deltaC);
+	public void setYield(LandCoverType forestType, int age, double woodProd) {
+		AgeFunction ageFunc = woodYields.get(forestType);
+		if (ageFunc == null) {
+			ageFunc = new AgeFunction();
+			woodYields.put(forestType, ageFunc);
+		}
+		ageFunc.setData(age, woodProd);
 	}
 	
 	public void setDefaultForMissingData() {
 		for (LandCoverType forestType : LandCoverType.getNaturalTypes()) {
-			for (int t = 0; t <= 100; t += 5) {
-				setForestGrowthFunction(forestType, t, 0);
+			for (int t = 0; t <= 300; t += 1) {
+				setYield(forestType, t, 0);
 			}
 		}
 	}
 	
-	public double getGrowth(LandCoverType forestType, int time) {
-		int lowerTime = (time/5) * 5;
-		int upperTime = lowerTime + 5;
-		Double lowerVal = forestGrowthFunction.get(forestType, lowerTime);
-		Double upperVal = forestGrowthFunction.get(forestType, upperTime);
-		double factor = ((double)(time - lowerTime)) / (upperTime - lowerTime);
-		Double d = Interpolator.interpolate(lowerVal, upperVal, factor);
-		return d;
-	}
-	
-	public double getCumulGrowth(LandCoverType forestType, int startTime, int endTime) {
-		double cumulGrowth = 0;
-		for (int i=startTime; i<=endTime; i++) {
-			cumulGrowth += getGrowth(forestType, i);
-		}
-		return cumulGrowth;
+	public double getYield(LandCoverType forestType, int age) {
+		return woodYields.get(forestType).getValue(age);
 	}
-	
-	public Map<Integer, Double> getGrowthFunction(LandCoverType forestType) {
-		return forestGrowthFunction.getInnerMap(forestType);
+
+	public AgeFunction getWoodYields(LandCoverType forestType) {
+		return woodYields.get(forestType);
 	}
 	
 	public int getMaxAge(LandCoverType forestType) {
-		return Collections.max(forestGrowthFunction.getInnerMap(forestType).keySet());
+		return woodYields.get(forestType).getDomainLimit();
 	}
 }
diff --git a/src/ac/ed/lurg/forestry/ForestManager.java b/src/ac/ed/lurg/forestry/ForestManager.java
index ea0f67c640cf085e37e47994041632c495b92779..ece6e6a5f6c4134fe341c9cbf574dd6a291687ce 100644
--- a/src/ac/ed/lurg/forestry/ForestManager.java
+++ b/src/ac/ed/lurg/forestry/ForestManager.java
@@ -7,200 +7,69 @@ import java.util.HashMap;
 import java.util.Map;
 
 import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.carbon.CarbonFluxItem;
+import ac.ed.lurg.landuse.LandUseTile;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.utils.LogWriter;
 
 public class ForestManager implements Serializable {
 	private static final long serialVersionUID = 3093028776147903435L;
 	
-	private Map<LandCoverType, ArrayList<ForestStandItem>> forest;
-	private Map<LandCoverType, Double> forestArea;
-	private Map<LandCoverType, Double> unprotectedForestArea;
-	private ForestGrowthItem growthFunction;
+	private ForestGrowthItem woodYields;
 	private int optimalRotation;
 	
-	public ForestManager(Map<LandCoverType, Double> landCoverAreas, Map<LandCoverType, Double> unprotectedAreas, 
-			double initialForestedNaturalArea, ForestGrowthItem growthFunction) {
-		this.forest = new HashMap<LandCoverType, ArrayList<ForestStandItem>>();
-		this.forestArea = new HashMap<LandCoverType, Double>();
-		this.unprotectedForestArea = new HashMap<LandCoverType, Double>();
-		for (LandCoverType forestType : LandCoverType.getNaturalTypes()) {
-			forestArea.put(forestType, landCoverAreas.get(forestType));
-			unprotectedForestArea.put(forestType, unprotectedAreas.get(forestType));
-		}
+	public ForestManager(ForestGrowthItem growthFunction) {
 		if (growthFunction == null) {
 			growthFunction = new ForestGrowthItem();
 			growthFunction.setDefaultForMissingData();
-//			LogWriter.printlnWarning("Null growthFunction");
-		}
-		this.growthFunction = growthFunction;
-		initialiseForests(landCoverAreas, unprotectedAreas, initialForestedNaturalArea);
-	}
-	
-	public void initialiseForests(Map<LandCoverType, Double> landCoverAreas, Map<LandCoverType, Double> unprotectedAreas, double initialForestedNaturalArea) {
-		
-		// Initial optimal rotation
-		calculateOptimalRotation(ModelConfig.INIT_WOOD_PRICE, ModelConfig.DISCOUNT_RATE);
-		
-		// Initial managed forest. Assume stands evenly distributed.
-		ArrayList<ForestStandItem> forestStands = new ArrayList<ForestStandItem>();
-		LandCoverType forestType = LandCoverType.TIMBER_FOREST;
-		double unprotArea = unprotectedAreas.get(forestType);
-		double protArea = landCoverAreas.get(forestType) - unprotArea;
-		double standArea = unprotArea / (optimalRotation + 1);
-		
-		if (unprotArea > 0) {
-			for (int age = 0; age <= optimalRotation; age++) {
-				double yield = growthFunction.getCumulGrowth(forestType, 0, age);
-				forestStands.add(new ForestStandItem(forestType, standArea, age, yield, false));
-			}
-		}
-		// Protected managed forest
-		if (protArea > 0) {
-			forestStands.add(new ForestStandItem(forestType, protArea, ModelConfig.FOREST_INIT_AGE, 
-					growthFunction.getCumulGrowth(forestType, 0, ModelConfig.FOREST_INIT_AGE), true));
-		}
-		
-		if (!forestStands.isEmpty())
-			forest.put(forestType, forestStands);
-		
-		// Natural forests
-		ArrayList<ForestStandItem> naturalStands= new ArrayList<ForestStandItem>();
-		
-		double protectedFract = (landCoverAreas.get(LandCoverType.NATURAL) - unprotectedAreas.get(LandCoverType.NATURAL)) / landCoverAreas.get(LandCoverType.NATURAL);
-		double natForestedArea = initialForestedNaturalArea;
-		double natForestedProtectedArea = natForestedArea * protectedFract;
-		double natOtherArea = landCoverAreas.get(LandCoverType.NATURAL) - initialForestedNaturalArea;
-		double natOtherProtectedArea = natOtherArea * protectedFract;
-		double yield = growthFunction.getCumulGrowth(LandCoverType.NATURAL, 0, ModelConfig.FOREST_INIT_AGE);
-		
-		if (natForestedArea - natForestedProtectedArea > 0)
-			naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natForestedArea - natForestedProtectedArea, ModelConfig.FOREST_INIT_AGE, yield, false));
-		
-		if (natForestedProtectedArea > 0)
-			naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natForestedProtectedArea, ModelConfig.FOREST_INIT_AGE, yield, true));
-		
-		if (natOtherArea - natOtherProtectedArea > 0)
-			naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natOtherArea - natOtherProtectedArea, 0, 0.0, false));
-		
-		if (natOtherProtectedArea > 0)
-			naturalStands.add(new ForestStandItem(LandCoverType.NATURAL, natOtherProtectedArea, 0, 0.0, true));	
-		
-		if (!naturalStands.isEmpty())
-			forest.put(LandCoverType.NATURAL, naturalStands);
-	}
-	
-	public void growForests() {
-		for (LandCoverType forestType : LandCoverType.getNaturalTypes()) {
-			ArrayList<ForestStandItem> forestStands = forest.get(forestType);
-			if (forestStands == null) 
-				continue;
-			for (ForestStandItem stand: forestStands) {
-				stand.updateAgeAndYield(growthFunction);
-			}
 		}
+		this.woodYields = growthFunction;
 	}
 	
 	public void setGrowthFunction(ForestGrowthItem growthFunction) {
-		this.growthFunction = growthFunction;
+		this.woodYields = growthFunction;
 	}
 	
 	public void calculateOptimalRotation(double woodPrice, double discountRate) {
 		Map<Double, Integer> levMap = new HashMap<Double, Integer>(); // Land Expected Value
-		int maxAge = growthFunction.getMaxAge(LandCoverType.TIMBER_FOREST);
+		int maxAge = woodYields.getMaxAge(LandCoverType.TIMBER_FOREST);
 		
 		for (int age = 0; age <= maxAge; age++) {
-			double yield = growthFunction.getCumulGrowth(LandCoverType.TIMBER_FOREST, 0, age);
+			double yield = woodYields.getYield(LandCoverType.TIMBER_FOREST, age);
 			double lev = (woodPrice * yield * Math.exp(-discountRate * age) - ModelConfig.FOREST_ESTABLISHMENT_COST) / (1 - Math.exp(-discountRate * age));
 			levMap.put(lev, age);
 		}
 		// TODO potential edge case with tied LEV
 		double maxLev = Collections.max(levMap.keySet());
 		optimalRotation = levMap.get(maxLev);
-		if (optimalRotation < 100) {
-			int foo = 1;
-		}
 	}
 	
 	public int getOptimalRotation() {
 		return optimalRotation;
 	}
 	
-	public double harvestWood() {
-		ArrayList<ForestStandItem> timberForest = forest.get(LandCoverType.TIMBER_FOREST);
-		if (timberForest == null)
-			return 0;
-		
-		double woodHarvested = 0;
-		for (ForestStandItem stand : timberForest) {
-			if (stand.getAge() >= optimalRotation) {
-				double standHarvest = stand.harvestStandWood();
-				woodHarvested += standHarvest;
-			}
-		}
-		return woodHarvested;
-	}
-	
-	public void moveForestArea(LandCoverType toType, LandCoverType fromType, double area) {
-		
-		ArrayList<ForestStandItem> toForest = forest.get(toType);
-		ArrayList<ForestStandItem> fromForest = forest.get(fromType);
-		
-		// Edge case where both types are not forest
-		if (toForest == null && fromForest == null) {
-			return;
-		} else
-		
-		// Moving from forest to non-forest i.e. removing forest
-		if (fromForest != null) {
-			double totalCurrentArea = unprotectedForestArea.get(fromType);
-			for (ForestStandItem s: fromForest) {
-				if (s.isProtected())
-					continue;
-				double areaToRemove = s.getArea() / totalCurrentArea * area; 
-				s.removeArea(areaToRemove);
-			}
-			unprotectedForestArea.put(fromType, unprotectedForestArea.get(fromType) - area);
-			forestArea.put(fromType, forestArea.get(fromType) - area);
-		}
-		
-		// Moving from non-forest to forest i.e. adding new forest
-		if (toForest != null) {
-			ForestStandItem newStand = new ForestStandItem(toType, area, 0, 0, false);
-			forest.get(toType).add(newStand);
-			unprotectedForestArea.put(toType, unprotectedForestArea.get(toType) + area);
-			forestArea.put(toType, forestArea.get(toType) + area);
-		}		
+	public double getStandYield(LandUseTile tile) {
+		return woodYields.getYield(tile.getLcType(), tile.getAge());
 	}
 	
-	public double getAverageYield(LandCoverType forestType) {
-		ArrayList<ForestStandItem> stands = forest.get(forestType);
-		if (stands == null) 
+	public double getAverageYield(ArrayList<LandUseTile> forestTiles) {
+		if (forestTiles == null) 
 			return 0;
-		Double totalWoodStock = stands.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getWoodStock()).sum();
+		Double totalWoodStock = forestTiles.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getArea() * getStandYield(o)).sum();
+		if (totalWoodStock > 0) {
+			int foo = 1;
+		}
+		double totalArea = forestTiles.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getArea()).sum();
 		if (totalWoodStock == null || totalWoodStock == 0) {
 			return 0.0;
 		} else {
-			return totalWoodStock / unprotectedForestArea.get(forestType);
+			return totalWoodStock / totalArea;
 		}	
 	}
 	
 	public double getYieldAtRotation() {
-		return growthFunction.getCumulGrowth(LandCoverType.TIMBER_FOREST, 0, optimalRotation);
+		return woodYields.getYield(LandCoverType.TIMBER_FOREST, optimalRotation);
 	}
 	
-	public void runAreaCheck() {
-		for (Map.Entry<LandCoverType, ArrayList<ForestStandItem>> entry : forest.entrySet()) {
-			LandCoverType forestType = entry.getKey();
-			ArrayList<ForestStandItem> stands = entry.getValue();
-			double totalArea = stands.stream().mapToDouble(o -> o.getArea()).sum();
-			double totalUnprotectedArea = stands.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getArea()).sum();
-			if (totalArea != forestArea.get(forestType)) {
-				LogWriter.printlnError("ForestManager: total forest area and stand area different: " + forestArea.get(forestType) + ", " + totalArea);
-			}
-			if (totalUnprotectedArea != unprotectedForestArea.get(forestType)) {
-				LogWriter.printlnError("ForestManager: total forest area and stand area different: " + unprotectedForestArea.get(forestType) + ", " + totalUnprotectedArea);
-			}
-		}
-	}
+
 }
diff --git a/src/ac/ed/lurg/forestry/ForestStandItem.java b/src/ac/ed/lurg/forestry/ForestStandItem.java
deleted file mode 100644
index 4fd962b09ca0be581b3ae274e5b491df62cb4461..0000000000000000000000000000000000000000
--- a/src/ac/ed/lurg/forestry/ForestStandItem.java
+++ /dev/null
@@ -1,80 +0,0 @@
-package ac.ed.lurg.forestry;
-
-import java.io.Serializable;
-
-import ac.ed.lurg.ModelConfig;
-import ac.ed.lurg.types.LandCoverType;
-import ac.ed.lurg.utils.LogWriter;
-
-public class ForestStandItem implements Serializable {
-	private static final long serialVersionUID = -5437305221519749016L;
-	
-	private LandCoverType forestType;
-	private double area;
-	private int age;
-	private double yield;
-	private boolean isProtected;
-	
-	public ForestStandItem(LandCoverType forestType, double area, int age, double yield, boolean isProtected) {
-		this.forestType = forestType;
-		this.area = area;
-		this.age = age;
-		this.yield = yield;
-		this.isProtected = isProtected;
-	}
-
-	public LandCoverType getForestType() {
-		return forestType;
-	}
-
-	public double getArea() {
-		return area;
-	}
-
-	public double getAge() {
-		return age;
-	}
-
-	public void setAge(int age) {
-		this.age = age;
-	}
-
-	public double getYield() {
-		return yield;
-	}
-	
-	public void setYield(double yield) {
-		this.yield = yield;
-	}
-
-	public boolean isProtected() {
-		return isProtected;
-	}
-
-	public void removeArea(double a) {
-		if (a > area) {
-			area -= Math.min(area, a);
-			LogWriter.printlnError("ForestItem: attempting to remove too much area");
-		} else {
-			area -= a;
-		}	
-	}
-	
-	public void updateAgeAndYield(ForestGrowthItem forestGrowthFunction) {
-		age += ModelConfig.TIMESTEP_SIZE;
-		yield += forestGrowthFunction.getGrowth(forestType, age);
-		
-	}
-	
-	public double harvestStandWood() {
-		double harvest = getArea() * getYield();
-		setAge(0);
-		setYield(0);
-		return harvest;
-	}
-	
-	public double getWoodStock() {
-		return area * yield;
-	}
-
-}
diff --git a/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java
index f2de87ce9e96c0457ef57ed64a3791ed669998bd..5ff50fe5322b988fce0affb3c7690c88a35c804e 100644
--- a/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java
+++ b/src/ac/ed/lurg/forestry/NaturalForestGrowthDataReader.java
@@ -16,27 +16,27 @@ public class NaturalForestGrowthDataReader extends AbstractTabularRasterReader<F
 	}
 
 	protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 0, 0.0);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 5,  getValueForCol(rowValues, "5") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 10, getValueForCol(rowValues, "10") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 15, getValueForCol(rowValues, "15") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 20, getValueForCol(rowValues, "20") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 25, getValueForCol(rowValues, "25") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 30, getValueForCol(rowValues, "30") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 35, getValueForCol(rowValues, "35") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 40, getValueForCol(rowValues, "40") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 45, getValueForCol(rowValues, "45") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 50, getValueForCol(rowValues, "50") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 55, getValueForCol(rowValues, "55") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 60, getValueForCol(rowValues, "60") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 65, getValueForCol(rowValues, "65") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 70, getValueForCol(rowValues, "70") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 75, getValueForCol(rowValues, "75") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 80, getValueForCol(rowValues, "80") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 85, getValueForCol(rowValues, "85") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 90, getValueForCol(rowValues, "90") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 95, getValueForCol(rowValues, "95") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.NATURAL, 100, getValueForCol(rowValues, "100") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 0, 0.0);
+		item.setYield(LandCoverType.NATURAL, 5,  getValueForCol(rowValues, "5") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 10, getValueForCol(rowValues, "10") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 15, getValueForCol(rowValues, "15") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 20, getValueForCol(rowValues, "20") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 25, getValueForCol(rowValues, "25") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 30, getValueForCol(rowValues, "30") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 35, getValueForCol(rowValues, "35") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 40, getValueForCol(rowValues, "40") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 45, getValueForCol(rowValues, "45") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 50, getValueForCol(rowValues, "50") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 55, getValueForCol(rowValues, "55") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 60, getValueForCol(rowValues, "60") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 65, getValueForCol(rowValues, "65") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 70, getValueForCol(rowValues, "70") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 75, getValueForCol(rowValues, "75") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 80, getValueForCol(rowValues, "80") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 85, getValueForCol(rowValues, "85") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 90, getValueForCol(rowValues, "90") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 95, getValueForCol(rowValues, "95") * conversionFactor);
+		item.setYield(LandCoverType.NATURAL, 100, getValueForCol(rowValues, "100") * conversionFactor);
 
                                                                         
 	}
diff --git a/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java b/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java
index f2d8135d474dd3290840ec9c1a92fa55254795ce..aaa8d2527f268232e310f7155c33cb0d70d1e2fb 100644
--- a/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java
+++ b/src/ac/ed/lurg/forestry/TimberForestGrowthDataReader.java
@@ -16,27 +16,27 @@ public class TimberForestGrowthDataReader extends AbstractTabularRasterReader<Fo
 	}
 
 	protected void setData(RasterKey key, ForestGrowthItem item, Map<String, Double> rowValues) {
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 0, 0.0);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 5,  getValueForCol(rowValues, "5") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
-		item.setForestGrowthFunction(LandCoverType.TIMBER_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 0, 0.0);
+		item.setYield(LandCoverType.TIMBER_FOREST, 5,  getValueForCol(rowValues, "5") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 10, getValueForCol(rowValues, "10") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 15, getValueForCol(rowValues, "15") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 20, getValueForCol(rowValues, "20") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 25, getValueForCol(rowValues, "25") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 30, getValueForCol(rowValues, "30") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 35, getValueForCol(rowValues, "35") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 40, getValueForCol(rowValues, "40") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 45, getValueForCol(rowValues, "45") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 50, getValueForCol(rowValues, "50") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 55, getValueForCol(rowValues, "55") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 60, getValueForCol(rowValues, "60") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 65, getValueForCol(rowValues, "65") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 70, getValueForCol(rowValues, "70") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 75, getValueForCol(rowValues, "75") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 80, getValueForCol(rowValues, "80") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 85, getValueForCol(rowValues, "85") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 90, getValueForCol(rowValues, "90") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 95, getValueForCol(rowValues, "95") * conversionFactor);
+		item.setYield(LandCoverType.TIMBER_FOREST, 100, getValueForCol(rowValues, "100") * conversionFactor);
 
                                                                         
 	}
diff --git a/src/ac/ed/lurg/landuse/CarbonFluxItem.java b/src/ac/ed/lurg/landuse/CarbonFluxItem.java
deleted file mode 100644
index ce7452d61ef1f98bdc2219cbe10043514e9e9eef..0000000000000000000000000000000000000000
--- a/src/ac/ed/lurg/landuse/CarbonFluxItem.java
+++ /dev/null
@@ -1,43 +0,0 @@
-package ac.ed.lurg.landuse;
-
-import ac.sac.raster.RasterItem;
-import java.util.HashMap;
-import java.util.Map;
-
-import ac.ed.lurg.types.LandCoverType;
-
-public class CarbonFluxItem implements RasterItem {
-	Map<LandCoverType, Map<LandCoverType, Double>> carbonFluxes = new HashMap<LandCoverType, Map<LandCoverType, Double>>(); 
-	
-	public void setCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover, double carbonFlux) {
-		
-		if (carbonFluxes.containsKey(previousLandCover)) {
-			carbonFluxes.get(previousLandCover).put(newLandCover, carbonFlux);
-		} else {
-			Map<LandCoverType, Double> cfMap = new HashMap<LandCoverType, Double>();
-			cfMap.put(newLandCover, carbonFlux);
-			carbonFluxes.put(previousLandCover, cfMap);
-		}
-	}
-	
-	public double getCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover) {
-		return carbonFluxes.get(previousLandCover).get(newLandCover);
-	}
-	
-	public boolean checkForKeys(LandCoverType previousLandCover, LandCoverType newLandCover) {
-		if (carbonFluxes.containsKey(previousLandCover)) {
-			if (carbonFluxes.get(previousLandCover).containsKey(newLandCover)) {
-				return true;
-			} else {
-				return false;
-			}
-		} else {
-			return false;
-		}
-	}
-	
-	public Map<LandCoverType, Map<LandCoverType, Double>> getCarbonFluxes() {
-		return carbonFluxes;
-	}
-
-}
diff --git a/src/ac/ed/lurg/landuse/ConversionCostReader.java b/src/ac/ed/lurg/landuse/ConversionCostReader.java
index 05aec14ece830af942a69ede81c8dfc8afc08951..c8a28a6a3ac9f284eaea2d3f1df9f5705983e99f 100644
--- a/src/ac/ed/lurg/landuse/ConversionCostReader.java
+++ b/src/ac/ed/lurg/landuse/ConversionCostReader.java
@@ -3,10 +3,11 @@ package ac.ed.lurg.landuse;
 import java.io.BufferedReader;
 import java.io.FileReader;
 import java.io.IOException;
+import java.util.HashMap;
+import java.util.Map;
 
 import ac.ed.lurg.ModelConfig;
 import ac.ed.lurg.types.LandCoverType;
-import ac.ed.lurg.utils.DoubleMap;
 import ac.ed.lurg.utils.LogWriter;
 
 public class ConversionCostReader {
@@ -15,12 +16,12 @@ public class ConversionCostReader {
 	private static final int TO_COL = 1;
 	private static final int COST_COL = 2;
 	
-	public DoubleMap<LandCoverType, LandCoverType, Double> read() {
-		DoubleMap<LandCoverType, LandCoverType, Double> conversionCosts = new DoubleMap<LandCoverType, LandCoverType, Double>();
+	public Map<LccKey, Double> read() {
+		Map<LccKey, Double> conversionCosts = new HashMap<LccKey, Double>();
 		String filename = ModelConfig.CONVERSION_COST_FILE;
 		try {
 			BufferedReader reader = new BufferedReader(new FileReader(filename)); 
-			String line, fromName, toName;
+			String line, toName, fromName;
 			Double cost;
 			reader.readLine(); // read header
 
@@ -34,10 +35,10 @@ public class ConversionCostReader {
 				toName = tokens[TO_COL];
 				cost = Double.parseDouble(tokens[COST_COL]);
 				
-				LandCoverType fromLC = LandCoverType.getForName(fromName);
-				LandCoverType toLC = LandCoverType.getForName(toName);
+				LandCoverType fromLc = LandCoverType.getForName(fromName);
+				LandCoverType toLc = LandCoverType.getForName(toName);
 				
-				conversionCosts.put(fromLC, toLC, cost);
+				conversionCosts.put(new LccKey(fromLc, toLc), cost);
 			} 
 			reader.close(); 
 		
@@ -49,5 +50,38 @@ public class ConversionCostReader {
 		
 		return conversionCosts;
 	}
+	
+	public Map<LccKey, Double> calcFromConfig() {
+		Map<LccKey, Double> conversionCosts = new HashMap<LccKey, Double>();
+		
+		conversionCosts.put(new LccKey(LandCoverType.CROPLAND, LandCoverType.PASTURE), ModelConfig.CROP_DECREASE_COST + ModelConfig.PASTURE_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CROPLAND, LandCoverType.NATURAL), ModelConfig.CROP_DECREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST), ModelConfig.CROP_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST), ModelConfig.CROP_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
+		
+		conversionCosts.put(new LccKey(LandCoverType.PASTURE, LandCoverType.CROPLAND), ModelConfig.PASTURE_DECREASE_COST + ModelConfig.CROP_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.PASTURE, LandCoverType.NATURAL), ModelConfig.PASTURE_DECREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.PASTURE, LandCoverType.TIMBER_FOREST), ModelConfig.PASTURE_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.PASTURE, LandCoverType.CARBON_FOREST), ModelConfig.PASTURE_DECREASE_COST + ModelConfig.MANAGED_FOREST_INCREASE_COST);
+		
+		conversionCosts.put(new LccKey(LandCoverType.NATURAL, LandCoverType.CROPLAND), ModelConfig.AGRI_EXPANSION_COST_BASE + ModelConfig.CROP_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
+		conversionCosts.put(new LccKey(LandCoverType.NATURAL, LandCoverType.PASTURE), ModelConfig.AGRI_EXPANSION_COST_BASE + ModelConfig.PASTURE_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
+		conversionCosts.put(new LccKey(LandCoverType.NATURAL, LandCoverType.TIMBER_FOREST), ModelConfig.MANAGED_FOREST_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
+		conversionCosts.put(new LccKey(LandCoverType.NATURAL, LandCoverType.CARBON_FOREST), ModelConfig.MANAGED_FOREST_INCREASE_COST + ModelConfig.INFRASTRUCTURE_EXPANSION_COST);
+		
+		conversionCosts.put(new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CROPLAND), ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.PASTURE), ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.PASTURE_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.NATURAL), ModelConfig.MANAGED_FOREST_DECREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST), 0.0);
+		conversionCosts.put(new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CARBON_FOREST), ModelConfig.FOREST_MANAGEMENT_COST);
+		
+		conversionCosts.put(new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.CROPLAND), ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.CROP_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.PASTURE), ModelConfig.MANAGED_FOREST_DECREASE_COST + ModelConfig.AGRI_EXPANSION_COST_BASE_MANAGED_FOREST + ModelConfig.PASTURE_INCREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.NATURAL), ModelConfig.MANAGED_FOREST_DECREASE_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.TIMBER_FOREST), ModelConfig.FOREST_MANAGEMENT_COST);
+		conversionCosts.put(new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.CARBON_FOREST), ModelConfig.FOREST_MANAGEMENT_COST);
+		
+		return conversionCosts;
+	}
 
 }
diff --git a/src/ac/ed/lurg/landuse/LandCoverChangeItem.java b/src/ac/ed/lurg/landuse/LandCoverChangeItem.java
deleted file mode 100644
index cdde2f6e8e10c07cd7fb32c030f30a1fd430848d..0000000000000000000000000000000000000000
--- a/src/ac/ed/lurg/landuse/LandCoverChangeItem.java
+++ /dev/null
@@ -1,46 +0,0 @@
-package ac.ed.lurg.landuse;
-
-import java.util.HashMap;
-import java.util.Map;
-
-import ac.ed.lurg.types.LandCoverType;
-
-public class LandCoverChangeItem {
-	private Map<Integer, Map<LandCoverType, Map<LandCoverType, Double>>> landCoverChange = new HashMap<Integer, Map<LandCoverType, Map<LandCoverType, Double>>>();
-	
-	public LandCoverChangeItem() {}
-	
-	public void setLandCoverChange(LandCoverType fromLC, LandCoverType toLC, Integer locId, double change) {
-		
-		if (landCoverChange.containsKey(locId)) {
-			Map<LandCoverType, Map<LandCoverType, Double>> locMap = landCoverChange.get(locId);
-			
-			if (locMap.containsKey(fromLC)) {
-				locMap.get(fromLC).put(toLC, change);
-			} 
-			else {
-				Map<LandCoverType, Double> toMap = new HashMap<LandCoverType, Double>();
-				toMap.put(toLC, change);
-				
-				landCoverChange.get(locId).put(fromLC, toMap);
-			}
-		}
-		else {
-			Map<LandCoverType, Double> toMap = new HashMap<LandCoverType, Double>();
-			toMap.put(toLC, change);
-			
-			Map<LandCoverType, Map<LandCoverType, Double>> fromMap = new HashMap<LandCoverType, Map<LandCoverType, Double>>();
-			fromMap.put(fromLC, toMap);
-			
-			landCoverChange.put(locId, fromMap);
-		}
-	}
-	
-	public double getLandCoverChange(LandCoverType fromLC, LandCoverType toLC, Integer locId) {
-		return landCoverChange.get(locId).get(fromLC).get(toLC);
-	}
-	
-	public Map<LandCoverType, Map<LandCoverType, Double>> getLandCoverChangeForLocation(Integer locId) {
-		return landCoverChange.get(locId);
-	}
-}
diff --git a/src/ac/ed/lurg/landuse/LandUseItem.java b/src/ac/ed/lurg/landuse/LandUseItem.java
index 76f22b53404313c0f0fd9c48e559d3ceb71385c9..7333d94605cc31f50d0f4933194a9c06fa13258b 100644
--- a/src/ac/ed/lurg/landuse/LandUseItem.java
+++ b/src/ac/ed/lurg/landuse/LandUseItem.java
@@ -1,6 +1,7 @@
 package ac.ed.lurg.landuse;
 
 import java.io.Serializable;
+import java.util.ArrayList;
 import java.util.Collection;
 import java.util.HashMap;
 import java.util.Map;
@@ -21,6 +22,7 @@ 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>(); //Mha
 	private Map<LandCoverType, Double> unprotectedAreas = new HashMap<LandCoverType, Double>(); //Mha
+	private Map<LandCoverType, ArrayList<LandUseTile>> landUseTiles = new HashMap<LandCoverType, ArrayList<LandUseTile>>();
 	private double protectedArea; //protected area in Mha
 	private double unavailableArea; //area unavailable due to altitude etc 
 	private double suitableArea; //Mha
@@ -38,7 +40,10 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 			setCropFraction(CropType.MAIZE, 0.5);
 			setUnavailableArea(landCover.getUnavailableFract());
 			updateSuitableArea(ModelConfig.BASE_YEAR);
-			forestManager = new ForestManager(landCoverAreas, unprotectedAreas, landCover.getInitialForestedNaturalArea(), growthFunction);
+			forestManager = new ForestManager(growthFunction);
+			initialiseForests(landCover.getInitialForestedNaturalArea());
+			initialiseAgriLand();
+			runAreaCheck();
 		}
 	}
 	
@@ -48,11 +53,129 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 		cropFractions.putAll(luItemToCopy.cropFractions);
 		landCoverAreas.putAll(luItemToCopy.landCoverAreas);
 		unprotectedAreas.putAll(luItemToCopy.unprotectedAreas);
+		landUseTiles.putAll(luItemToCopy.landUseTiles);
 		protectedArea = (luItemToCopy.protectedArea);
 		suitableArea = (luItemToCopy.suitableArea);
 		forestManager = (luItemToCopy.forestManager);
 	}
 	
+	private void initialiseForests(double initialForestedNaturalArea) {
+		
+		// Initial optimal rotation
+		forestManager.calculateOptimalRotation(ModelConfig.INIT_WOOD_PRICE, ModelConfig.DISCOUNT_RATE);
+		int optimalRotation = forestManager.getOptimalRotation();
+		
+		// Initial managed forest. Assume stands evenly distributed up to rotation age.
+		ArrayList<LandUseTile> forestStands = new ArrayList<LandUseTile>();
+		LandCoverType forestType = LandCoverType.TIMBER_FOREST;
+		double unprotArea = unprotectedAreas.get(forestType);
+		double protArea = landCoverAreas.get(forestType) - unprotArea;
+		double standArea = unprotArea / (optimalRotation + 1);
+		
+		if (unprotArea > 0) {
+			for (int age = 0; age <= optimalRotation; age++) {
+				forestStands.add(new LandUseTile(forestType, standArea, age, false));
+			}
+		}
+		// Protected managed forest
+		if (protArea > 0) {
+			forestStands.add(new LandUseTile(forestType, protArea, ModelConfig.FOREST_INIT_AGE, true));
+		}
+		
+		landUseTiles.put(forestType, forestStands);
+		landUseTiles.put(LandCoverType.CARBON_FOREST, new ArrayList<LandUseTile>());
+		
+		// Natural forests
+		ArrayList<LandUseTile> naturalStands= new ArrayList<LandUseTile>();
+		
+		double protectedFract = (landCoverAreas.get(LandCoverType.NATURAL) - unprotectedAreas.get(LandCoverType.NATURAL)) / landCoverAreas.get(LandCoverType.NATURAL);
+		double natForestedArea = initialForestedNaturalArea;
+		double natForestedProtectedArea = natForestedArea * protectedFract;
+		double natOtherArea = landCoverAreas.get(LandCoverType.NATURAL) - initialForestedNaturalArea;
+		double natOtherProtectedArea = natOtherArea * protectedFract;
+		
+		if (natForestedArea - natForestedProtectedArea > 0)
+			naturalStands.add(new LandUseTile(LandCoverType.NATURAL, natForestedArea - natForestedProtectedArea, ModelConfig.FOREST_INIT_AGE, false));
+		
+		if (natForestedProtectedArea > 0)
+			naturalStands.add(new LandUseTile(LandCoverType.NATURAL, natForestedProtectedArea, ModelConfig.FOREST_INIT_AGE, true));
+		
+		if (natOtherArea - natOtherProtectedArea > 0)
+			naturalStands.add(new LandUseTile(LandCoverType.NATURAL, natOtherArea - natOtherProtectedArea, 0, false));
+		
+		if (natOtherProtectedArea > 0)
+			naturalStands.add(new LandUseTile(LandCoverType.NATURAL, natOtherProtectedArea, 0, true));	
+		
+		landUseTiles.put(LandCoverType.NATURAL, naturalStands);
+	}
+	
+	private void initialiseAgriLand() {
+		ArrayList<LandUseTile> cropTiles = new ArrayList<LandUseTile>();
+		cropTiles.add(new LandUseTile(LandCoverType.CROPLAND, unprotectedAreas.get(LandCoverType.CROPLAND), 30, false));
+		landUseTiles.put(LandCoverType.CROPLAND, cropTiles);
+		
+		ArrayList<LandUseTile> pastTiles = new ArrayList<LandUseTile>();
+		pastTiles.add(new LandUseTile(LandCoverType.PASTURE, unprotectedAreas.get(LandCoverType.PASTURE), 30, false));
+		landUseTiles.put(LandCoverType.PASTURE, pastTiles);
+	}
+	
+	public void moveTilesArea(LandCoverType toType, LandCoverType fromType, double area) {
+		// Note if toType == fromType, area will be removed from fromType and a new tile will be created with age 0
+		
+		if (area == 0)
+			return;
+		
+		ArrayList<LandUseTile> toTiles = landUseTiles.get(toType);
+		ArrayList<LandUseTile> fromTiles = landUseTiles.get(fromType);
+		
+		if (fromTiles == null) {
+			LogWriter.printlnError("LandUseItem: cannot move area between tiles as source is null");
+			return;
+		}
+	
+		// Remove area from tiles
+			double totalCurrentArea = fromTiles.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getArea()).sum();
+			for (LandUseTile s: fromTiles) {
+				if (s.isProtected())
+					continue;
+
+				double areaToRemove = s.getArea() / totalCurrentArea * area; 
+				s.removeArea(areaToRemove);
+			}		
+		
+		// Add new tiles	
+			LandUseTile newStand = new LandUseTile(toType, area, 0, false);
+			toTiles.add(newStand);	
+			
+			runAreaCheck();
+	}
+	
+	public void runAreaCheck() {
+		double threshold = 1e-10;
+		for (Map.Entry<LandCoverType, ArrayList<LandUseTile>> entry : landUseTiles.entrySet()) {
+			LandCoverType lcType = entry.getKey();
+			ArrayList<LandUseTile> tiles = entry.getValue();
+			double totalArea = tiles.stream().mapToDouble(o -> o.getArea()).sum();
+			double totalUnprotectedArea = tiles.stream().filter(o -> !o.isProtected()).mapToDouble(o -> o.getArea()).sum();
+			if (totalArea - landCoverAreas.get(lcType) > threshold) {
+				LogWriter.printlnError("LandUseItem: total area and tiles area different: " + lcType.getName() + " " + landCoverAreas.get(lcType) + ", " + totalArea);
+			}
+			if (totalUnprotectedArea - unprotectedAreas.get(lcType) > threshold) {
+				LogWriter.printlnError("LandUseItem: unprotected area and tiles area different: " + lcType.getName() + " " + unprotectedAreas.get(lcType) + ", " + totalUnprotectedArea);
+			}
+		}
+	}
+	
+	public void updateTilesAge() {
+		for (Map.Entry<LandCoverType, ArrayList<LandUseTile>> entry : landUseTiles.entrySet()) {
+			LandCoverType lcType = entry.getKey();
+			ArrayList<LandUseTile> tiles = entry.getValue();
+			for (LandUseTile t : tiles) {
+				t.incrementAge();
+			}
+		}
+	}
+	
 	public void setProtectedFract(double protFrac) {
 		if (!Double.isNaN(ModelConfig.CONSTANT_PROTECTED_AREA_RATE))
 			protFrac=ModelConfig.CONSTANT_PROTECTED_AREA_RATE;
@@ -237,9 +360,8 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
         }
         */
         
-        if (toType.isNatural() || fromType.isNatural()) {
-        	forestManager.moveForestArea(toType, fromType, areaMoved);
-        }
+        moveTilesArea(toType, fromType, areaMoved);
+        
         
         return area - areaMoved;
 	}
@@ -392,15 +514,6 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 			suitableArea = currentAgri + getTotalUnprotectedNatural();  
 	}
 	
-	/*
-	public double getForestManagedFraction() {
-		double managed = getLandCoverArea(LandCoverType.TIMBER_FOREST) + getLandCoverArea(LandCoverType.CARBON_FOREST);
-		double unmanaged = getLandCoverArea(LandCoverType.UNMANAGED_FOREST);
-		double d = managed / (managed + unmanaged);
-		return (Double.isNaN(d) || Double.isInfinite(d)) ? 1.0 : d;
-	}
-	*/
-	
 	public CropToDouble getCropChanges(LandUseItem prevAreaAggItem) {
 		CropToDouble changes = new CropToDouble();
 		
@@ -416,16 +529,12 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 		forestManager.setGrowthFunction(newGrowthFunction);
 	}
 	
-	public void growForests() {
-		forestManager.growForests();
-	}
-	
 	public void calculateOptimalRotation(double woodPrice, double discountRate) {
 		forestManager.calculateOptimalRotation(woodPrice, discountRate);
 	}
 	
 	public double getWoodYield(LandCoverType forestType) {
-		return forestManager.getAverageYield(forestType);
+		return forestManager.getAverageYield(landUseTiles.get(forestType));
 	}
 	
 	public double getWoodYieldAtRotation() {
@@ -433,7 +542,19 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 	}
 	
 	public double getWoodHarvestFromRotation() {
-		return forestManager.harvestWood();
+		ArrayList<LandUseTile> tiles = landUseTiles.get(LandCoverType.TIMBER_FOREST);
+		if (tiles == null) 
+			return 0;
+		double totalHarvest = 0;
+		int rotationAge = forestManager.getOptimalRotation();
+		for (LandUseTile t : tiles) {
+			if (t.getAge() < rotationAge)
+				continue;
+			double yield = forestManager.getStandYield(t);
+			totalHarvest += yield * t.getArea();
+			t.resetAge(); // Stand age to 0
+		}
+		return totalHarvest;
 	}
 	
 	public int getForestRotation() {
diff --git a/src/ac/ed/lurg/landuse/LandUseTile.java b/src/ac/ed/lurg/landuse/LandUseTile.java
new file mode 100644
index 0000000000000000000000000000000000000000..9c826e6d5e48a18217a574d0d5119cb8281e0ff4
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/LandUseTile.java
@@ -0,0 +1,56 @@
+package ac.ed.lurg.landuse;
+
+import java.io.Serializable;
+
+import ac.ed.lurg.ModelConfig;
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.LogWriter;
+
+public class LandUseTile implements Serializable {
+	private static final long serialVersionUID = 3208516181615267472L;
+	
+	private LandCoverType lcType;
+	private double area;
+	private int age;
+	private boolean isProtected;
+	
+	public LandUseTile(LandCoverType lcType, double area, int age, boolean isProtected) {
+		this.lcType = lcType;
+		this.area = area;
+		this.age = age;
+		this.isProtected = isProtected;
+	}
+
+	public LandCoverType getLcType() {
+		return lcType;
+	}
+
+	public double getArea() {
+		return area;
+	}
+
+	public int getAge() {
+		return age;
+	}
+	
+	public boolean isProtected() {
+		return isProtected;
+	}
+
+	public void incrementAge() {
+		age += ModelConfig.TIMESTEP_SIZE;
+	}
+	
+	public void resetAge() {
+		age = 0;
+	}
+
+	public void removeArea(double a) {
+		if (a > area) {
+			LogWriter.printlnError("LandUseItem: attempting to remove too much area: " + a + ", available: " + area);
+			area -= Math.min(area, a);
+		} else {
+			area -= a;
+		}	
+	}
+}
diff --git a/src/ac/ed/lurg/landuse/LccKey.java b/src/ac/ed/lurg/landuse/LccKey.java
new file mode 100644
index 0000000000000000000000000000000000000000..5afbadb9443a1986261ea25a65f130950a2d39cf
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/LccKey.java
@@ -0,0 +1,47 @@
+package ac.ed.lurg.landuse;
+
+import ac.ed.lurg.types.LandCoverType;
+
+public class LccKey {
+	// Land Cover Change Key
+	private LandCoverType fromLc;
+	private LandCoverType toLc;
+	
+	public LccKey(LandCoverType fromLc, LandCoverType toLc) {
+		this.fromLc = fromLc;
+		this.toLc = toLc;
+	}
+
+	public LandCoverType getFromLc() {
+		return fromLc;
+	}
+
+	public LandCoverType getToLc() {
+		return toLc;
+	}
+
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		result = prime * result + ((fromLc == null) ? 0 : fromLc.hashCode());
+		result = prime * result + ((toLc == null) ? 0 : toLc.hashCode());
+		return result;
+	}
+
+	@Override
+	public boolean equals(Object obj) {
+		if (this == obj)
+			return true;
+		if (obj == null)
+			return false;
+		if (getClass() != obj.getClass())
+			return false;
+		LccKey other = (LccKey) obj;
+		if (fromLc != other.fromLc)
+			return false;
+		if (toLc != other.toLc)
+			return false;
+		return true;
+	}	
+}
diff --git a/src/ac/ed/lurg/utils/DualKey.java b/src/ac/ed/lurg/utils/DualKey.java
new file mode 100644
index 0000000000000000000000000000000000000000..0841eaedc6dfedaf090e18b673c16f2100e5849a
--- /dev/null
+++ b/src/ac/ed/lurg/utils/DualKey.java
@@ -0,0 +1,43 @@
+package ac.ed.lurg.utils;
+
+public class DualKey<K, Q> {
+	private K key1;
+	private Q key2;
+	
+	public DualKey(K key1, Q key2) {
+		this.key1 = key1;
+		this.key2 = key2;
+	}
+
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		result = prime * result + ((key1 == null) ? 0 : key1.hashCode());
+		result = prime * result + ((key2 == null) ? 0 : key2.hashCode());
+		return result;
+	}
+
+	@SuppressWarnings("unchecked")
+	@Override
+	public boolean equals(Object obj) {
+		if (this == obj)
+			return true;
+		if (obj == null)
+			return false;
+		if (getClass() != obj.getClass())
+			return false;
+		DualKey<K, Q> other = (DualKey<K, Q>) obj;
+		if (key1 == null) {
+			if (other.key1 != null)
+				return false;
+		} else if (!key1.equals(other.key1))
+			return false;
+		if (key2 == null) {
+			if (other.key2 != null)
+				return false;
+		} else if (!key2.equals(other.key2))
+			return false;
+		return true;
+	}	
+}