From e883162ff644b236adf816fedb446ee3759d2d04 Mon Sep 17 00:00:00 2001
From: Peter Alexander <peter@blackhillock.co.uk>
Date: Wed, 22 Jul 2015 10:59:18 +0100
Subject: [PATCH] Change values held by AreasItem
---
GAMS/IntExtOpt.gms | 3 +
src/ac/ed/lurg/ModelMain.java | 4 +-
src/ac/ed/lurg/country/CountryAgent.java | 15 +-
.../country/gams/GamsLocationOptimiser.java | 72 ++++---
.../lurg/country/gams/GamsLocationTest.java | 25 +--
.../country/gams/GamsRasterOptimiser.java | 192 ++++++++++--------
src/ac/ed/lurg/landuse/AreasItem.java | 77 ++++---
src/ac/ed/lurg/landuse/LandCoverItem.java | 8 +-
src/ac/ed/lurg/landuse/LandCoverReader.java | 6 +-
src/ac/ed/lurg/types/CropToDouble.java | 9 -
src/ac/ed/lurg/types/CropType.java | 2 +-
.../{LandDataType.java => LandCoverType.java} | 5 +-
src/ac/sac/raster/RasterSet.java | 1 +
13 files changed, 238 insertions(+), 181 deletions(-)
rename src/ac/ed/lurg/types/{LandDataType.java => LandCoverType.java} (79%)
diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index 5cef5c68..de004b9d 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -215,6 +215,9 @@ $gdxin
totalProdCost('meat') = sum(feed_crop, feedEnergy(feed_crop))
+ parameter totalCropland(location);
+ totalCropland(location) = sum(crop_less_pasture, area.l(crop_less_pasture, location));
+
display totalProdCost, totalProd, totalImportEnergy, feedEnergy, totalProdCost
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 711672f0..9500e56e 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -28,7 +28,7 @@ import ac.ed.lurg.output.RasterOutputer;
import ac.ed.lurg.types.CropToDoubleMap;
import ac.ed.lurg.types.CropType;
import ac.ed.lurg.types.FertiliserRate;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.ed.lurg.types.ModelFitType;
import ac.ed.lurg.utils.LogWriter;
import ac.ed.lurg.yield.LPJYieldResponseMapReader;
@@ -281,7 +281,7 @@ public class ModelMain {
String rootDir = ModelConfig.INITAL_LAND_COVER_DIR + File.separator;
- for (LandDataType landType : LandDataType.values()) {
+ for (LandCoverType landType : LandCoverType.values()) {
LandCoverReader lcReader = new LandCoverReader(initLC, landType);
initLC = lcReader.getRasterDataFromFile(rootDir + landType.getFileName());
}
diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java
index d23dc87d..4b2507f1 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -15,7 +15,6 @@ import ac.ed.lurg.landuse.IrrigationCostItem;
import ac.ed.lurg.landuse.LandCoverItem;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
-import ac.ed.lurg.types.LandDataType;
import ac.ed.lurg.utils.LogWriter;
import ac.ed.lurg.yield.YieldRaster;
import ac.sac.raster.RasterKey;
@@ -48,21 +47,15 @@ public class CountryAgent {
private RasterSet<AreasItem> convertInitialLC(RasterSet<LandCoverItem> initialLC) {
RasterSet<AreasItem> cropAreaRaster = new RasterSet<AreasItem>(initialLC.getHeaderDetails());
-
for (Map.Entry<RasterKey, LandCoverItem> entry : initialLC.entrySet()) {
AreasItem areasItem = new AreasItem();
RasterKey key = entry.getKey();
- LandCoverItem landCover = entry.getValue();
-
- if (landCover != null) {
- areasItem.setCropArea(CropType.WHEAT, landCover.getLandCover(LandDataType.CROPLAND)); // we allow free substitution between crop types so this is ok-ish
- areasItem.setCropArea(CropType.PASTURE, landCover.getLandCover(LandDataType.PASTURE));
- areasItem.setLandCoverArea(LandDataType.FOREST, landCover.getLandCover(LandDataType.FOREST));
- areasItem.setLandCoverArea(LandDataType.OTHER_NATURAL, landCover.getLandCover(LandDataType.OTHER_NATURAL));
- }
-
+ areasItem.setLandCoverAreas(entry.getValue());
+ areasItem.setCropFraction(CropType.WHEAT, 0.5); // random start, better if we could get data, but free substitution between crops so not critical
+ areasItem.setCropFraction(CropType.MAIZE, 0.5);
cropAreaRaster.put(key, areasItem);
}
+
return cropAreaRaster;
}
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index c04453ef..8bdf7f14 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -14,6 +14,8 @@ import ac.ed.lurg.landuse.Intensity;
import ac.ed.lurg.landuse.IrrigationCostItem;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.LazyHashMap;
import ac.ed.lurg.utils.LogWriter;
import ac.ed.lurg.yield.YieldResponsesItem;
@@ -92,21 +94,22 @@ public class GamsLocationOptimiser {
Integer locationId = entry.getKey();
String locString = Integer.toString(locationId);
AreasItem areasItem = entry.getValue();
+
+ double suitableLand = areasItem.getSuitableLand();
+ if (DEBUG) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, "suitableLand", suitableLand));
+ landP.addRecord(locString).setValue(suitableLand);
+
+ double pastureArea = areasItem.getLandCoverArea(LandCoverType.PASTURE);
+ if (DEBUG) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, "pastureArea", pastureArea));
+ setPreviousArea(prevCropP, locString, CropType.PASTURE.getGamsName(), pastureArea);
+ totalAgriLand += pastureArea;
- for (CropType cropType : CropType.values()) {
+ for (CropType cropType : CropType.getCropsLessPasture()) {
double d = areasItem.getCropArea(cropType);
if (DEBUG && d != 0) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, cropType.getGamsName(), d));
-
- Vector<String> v = new Vector<String>();
- v.add(cropType.getGamsName());
- v.add(locString);
- prevCropP.addRecord(v).setValue(d);
+ setPreviousArea(prevCropP, locString, cropType.getGamsName(), d);
totalAgriLand += d;
}
-
- double suitableLand = areasItem.getSuitableLand();
- if (DEBUG) LogWriter.println(String.format(" %d %15s,\t %.1f", locationId, "suitableLand", suitableLand));
- landP.addRecord(locString).setValue(suitableLand);
}
if (DEBUG) LogWriter.println(String.format(" Total agricultural\t %.1f", totalAgriLand));
@@ -183,6 +186,14 @@ public class GamsLocationOptimiser {
addItemMapParm(inDB.addParameter("cropAdj", 1), countryInput.getCropAdjustments());
}
+ private void setPreviousArea(GAMSParameter prevCropP, String locString, String cropTypeString, double d) {
+ Vector<String> v = new Vector<String>();
+ v.add(cropTypeString);
+ v.add(locString);
+ prevCropP.addRecord(v).setValue(d);
+ }
+
+ @SuppressWarnings("serial")
private GamsLocationOutput handleResults(GAMSDatabase outDB) {
int modelStatus = (int) outDB.getParameter("ms").findRecord().getValue();
System.out.println(
@@ -200,15 +211,20 @@ public class GamsLocationOptimiser {
GAMSParameter parmCropAdj = outDB.getParameter("cropAdj");
GAMSParameter parmProd = outDB.getParameter("totalProd");
GAMSParameter parmProdCost = outDB.getParameter("totalProdCost");
+ GAMSParameter parmCroplandArea = outDB.getParameter("totalCropland");
double totalArea = 0;
double area, fertIntensity, irrigIntensity, otherIntensity = Double.NaN, feedAmount, netImport, yield, unitEnergy, cropAdj, prod, prodCost;
- Map<Integer, IntensitiesItem> intensities = new HashMap<Integer, IntensitiesItem>();
- Map<Integer, AreasItem> cropAreas = new HashMap<Integer, AreasItem>();
+ LazyHashMap<Integer, IntensitiesItem> intensities = new LazyHashMap<Integer, IntensitiesItem>() {
+ protected IntensitiesItem createValue() { return new IntensitiesItem(); }
+ };
+ LazyHashMap<Integer, AreasItem> cropAreas = new LazyHashMap<Integer, AreasItem>() {
+ protected AreasItem createValue() { return new AreasItem(); }
+ };
+
Map<CropType, CropUsageData> cropUsageData = new HashMap<CropType, CropUsageData>();
Map<CropType, Double> cropAdjs = new HashMap<CropType, Double>();
- CropType prevCropType = null;
for (GAMSVariableRecord rec : varAreas) {
String itemName = rec.getKeys()[0];
@@ -223,39 +239,33 @@ public class GamsLocationOptimiser {
int locId = Integer.parseInt(locationName);
CropType cropType = CropType.getForGamsName(itemName);
- if (!cropType.equals(prevCropType)) {
+ if (!cropUsageData.containsKey(cropType)) { // then we must not have seen this crop type before, so need to do all non location specific stuff
feedAmount = varFeedAmount.findRecord(itemName).getLevel();
netImport = cropType.isImportedCrop() ? varNetImports.findRecord(itemName).getLevel() : 0;
cropAdj = getParmValue(parmCropAdj, itemName);
prod = getParmValue(parmProd, itemName);
prodCost = getParmValue(parmProdCost, itemName);
-
+
cropUsageData.put(cropType, new CropUsageData(feedAmount, netImport, prod, prodCost));
cropAdjs.put(cropType, cropAdj);
if (DEBUG) LogWriter.println(String.format("\n%s:\tfeedAmount= %.1f,\tnetImports= %.3f,\tprod= %.3f,\tprodCost= %.3f,\tcropAdj= %.3f", itemName, feedAmount, netImport, prod, prodCost, cropAdj));
}
if (area > 0) {
- if (DEBUG) LogWriter.println(String.format("\t location %s, %s:\tarea= %.1f,\tfert= %.3f,\tirrg= %.3f,\tintensity= %.3f",
- locationName, itemName, area, fertIntensity, irrigIntensity, otherIntensity));
-
- IntensitiesItem intensityItem = intensities.get(locId);
- if (intensityItem == null) {
- intensityItem = new IntensitiesItem();
- intensities.put(locId, intensityItem);
- }
+ if (DEBUG) LogWriter.println(String.format("\t location %s, %s:\tarea= %.1f,\tfert= %.3f,\tirrg= %.3f,\tintensity= %.3f", locationName, itemName, area, fertIntensity, irrigIntensity, otherIntensity));
+ IntensitiesItem intensityItem = intensities.lazyGet(locId);
intensityItem.setIntensity(cropType, new Intensity(fertIntensity, irrigIntensity, otherIntensity, yield, unitEnergy));
}
- AreasItem areaItem = cropAreas.get(locId);
- if (areaItem == null) {
- areaItem = new AreasItem();
- cropAreas.put(locId, areaItem);
- }
- areaItem.setCropArea(cropType, area);
- totalArea += area;
+ double croplandArea = getParmValue(parmCroplandArea, locationName);
+ AreasItem areasItem = cropAreas.lazyGet(locId);
+ if (cropType.equals(CropType.PASTURE))
+ areasItem.setLandCoverArea(LandCoverType.PASTURE, area);
+ else
+ areasItem.setLandCoverArea(LandCoverType.CROPLAND, croplandArea); // will set this multiple times, once for each arable crop, but doesn't really matter
- prevCropType = cropType;
+ areasItem.setCropFraction(cropType, croplandArea > 0 ? area/croplandArea : 0);
+ totalArea += area;
}
netImport = varNetImports.findRecord(CropType.MEAT.getGamsName()).getLevel();
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationTest.java b/src/ac/ed/lurg/country/gams/GamsLocationTest.java
index 73165127..84ec9637 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationTest.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationTest.java
@@ -8,7 +8,7 @@ import ac.ed.lurg.landuse.AreasItem;
import ac.ed.lurg.landuse.IrrigationCostItem;
import ac.ed.lurg.types.CommodityType;
import ac.ed.lurg.types.CropType;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.ed.lurg.types.YieldType;
import ac.ed.lurg.yield.YieldResponsesItem;
@@ -94,18 +94,19 @@ public class GamsLocationTest {
AreasItem getAreaItem(double cereal, double pasture) {
AreasItem dummyMap = new AreasItem();
- dummyMap.setCropArea(CropType.WHEAT, 5.0 + cereal);
- dummyMap.setCropArea(CropType.MAIZE, 5.0);
- dummyMap.setCropArea(CropType.RICE, 5.0);
- dummyMap.setCropArea(CropType.TROPICAL_CEREALS, 5.0);
- dummyMap.setCropArea(CropType.SOYBEAN, 5.0);
- dummyMap.setCropArea(CropType.OILCROPS, 5.0);
- dummyMap.setCropArea(CropType.PULSES, 5.0);
- dummyMap.setCropArea(CropType.STARCHY_ROOTS, 5.0);
- dummyMap.setCropArea(CropType.PASTURE, 5.0 + pasture);
+ dummyMap.setCropFraction(CropType.WHEAT, 0.2);
+ dummyMap.setCropFraction(CropType.MAIZE, 0.2);
+ dummyMap.setCropFraction(CropType.RICE, 0.1);
+ dummyMap.setCropFraction(CropType.TROPICAL_CEREALS, 0.1);
+ dummyMap.setCropFraction(CropType.SOYBEAN, 0.1);
+ dummyMap.setCropFraction(CropType.OILCROPS, 0.1);
+ dummyMap.setCropFraction(CropType.PULSES, 0.1);
+ dummyMap.setCropFraction(CropType.STARCHY_ROOTS, 0.1);
- dummyMap.setLandCoverArea(LandDataType.FOREST, 200 - pasture - cereal);
- dummyMap.setLandCoverArea(LandDataType.OTHER_NATURAL, 30);
+ dummyMap.setLandCoverArea(LandCoverType.PASTURE, pasture);
+ dummyMap.setLandCoverArea(LandCoverType.CROPLAND, cereal);
+ dummyMap.setLandCoverArea(LandCoverType.FOREST, 200 - pasture - cereal);
+ dummyMap.setLandCoverArea(LandCoverType.OTHER_NATURAL, 30);
return dummyMap;
}
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
index dcad1922..fac3693e 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
@@ -1,7 +1,6 @@
package ac.ed.lurg.country.gams;
import java.util.ArrayList;
-import java.util.Arrays;
import java.util.Collections;
import java.util.HashSet;
import java.util.List;
@@ -9,12 +8,13 @@ import java.util.Map;
import java.util.Map.Entry;
import java.util.Set;
+import ac.ed.lurg.ModelConfig;
import ac.ed.lurg.landuse.AreasItem;
import ac.ed.lurg.landuse.IntensitiesItem;
import ac.ed.lurg.landuse.IrrigationCostItem;
-import ac.ed.lurg.types.CropToDouble;
+import ac.ed.lurg.output.RasterOutputer;
import ac.ed.lurg.types.CropType;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.ed.lurg.types.YieldType;
import ac.ed.lurg.utils.LazyHashMap;
import ac.ed.lurg.utils.LogWriter;
@@ -61,9 +61,35 @@ public class GamsRasterOptimiser {
return new GamsRasterOutput(gamsOutput.getStatus(), newIntensityRaster, newAreaRaster, gamsOutput.getCommoditiesData(), gamsOutput.getCropAdjustments(), locationIdRaster);
}
+ private RasterSet<AreasItem> createWithSameLandCovers(RasterSet<AreasItem> toCopy) {
+ RasterSet<AreasItem> theCopy = new RasterSet<AreasItem>(toCopy.getHeaderDetails());
+
+ for (Entry<RasterKey, AreasItem> entry : toCopy.entrySet()) {
+ AreasItem newAreasItem = new AreasItem();
+ newAreasItem.setLandCoverAreas(entry.getValue());
+ theCopy.put(entry.getKey(), newAreasItem);
+ }
+
+ return theCopy;
+ }
+
+ @SuppressWarnings("unused")
+ private void checkedTotalAreas(RasterSet<AreasItem> areaRaster, String comment) {
+ for (LandCoverType l : LandCoverType.values()) {
+ double total = 0;
+ for (AreasItem a : areaRaster.values()) {
+ total += a.getLandCoverArea(l);
+ }
+
+ LogWriter.printlnError("Total Area " + comment + ": " + l.getName() + ": " + total);
+ }
+ }
+
private RasterSet<AreasItem> allocAreas(Map<Integer, ? extends AreasItem> prevAreasAgg, GamsLocationOutput gamsOutput) {
- RasterSet<AreasItem> prevAreaRaster = rasterInputData.getPreviousAreas();
- RasterSet<AreasItem> newAreaRaster = new RasterSet<AreasItem>(prevAreaRaster.getHeaderDetails());
+ RasterSet<AreasItem> newAreaRaster = createWithSameLandCovers(rasterInputData.getPreviousAreas());
+
+ //checkedTotalAreas(rasterInputData.getPreviousAreas(), "old");
+ //checkedTotalAreas(newAreaRaster, "new");
for (Map.Entry<Integer, AreasItem> entry : gamsOutput.getCropAreas().entrySet()) {
Integer locId = entry.getKey();
@@ -71,11 +97,13 @@ public class GamsRasterOptimiser {
AreasItem prevAreaAggItem = prevAreasAgg.get(locId);
Set<RasterKey> keys = mapping.get(locId);
- CropToDouble cropChangeTotals = newAreaAggItem.getCropChanges(prevAreaAggItem);
+ //checkedTotalAreas(newAreaRaster.popSubsetForKeys(new RasterSet<AreasItem>(), keys), locId + " setset");
+
+
// if (DEBUG) LogWriter.println("Processing location id " + locId);
- double pastureChange = cropChangeTotals.get(CropType.PASTURE);
- double croplandChange = cropChangeTotals.getCroplandTotal();
+ double pastureChange = newAreaAggItem.getLandCoverArea(LandCoverType.PASTURE) - prevAreaAggItem.getLandCoverArea(LandCoverType.PASTURE);
+ double croplandChange = newAreaAggItem.getLandCoverArea(LandCoverType.CROPLAND) - prevAreaAggItem.getLandCoverArea(LandCoverType.CROPLAND);
double prevForestToNaturalFraction = prevAreaAggItem.getForestToNaturalFraction();
double pastureFromCrop = 0;
@@ -84,103 +112,74 @@ public class GamsRasterOptimiser {
if (pastureChange > 0 && croplandChange < 0) {
pastureFromCrop = Math.min(pastureChange, -croplandChange);
- pastureFromNatural = (pastureChange > -croplandChange) ? 0 : pastureChange+croplandChange;
+ pastureFromNatural = (pastureChange > -croplandChange) ? pastureChange+croplandChange : 0;
+ cropFromNatural = croplandChange + pastureFromCrop;
}
else if (pastureChange < 0 && croplandChange > 0) {
pastureFromCrop = -Math.min(-pastureChange, croplandChange);
- cropFromNatural = (croplandChange > -pastureChange) ? 0 : -(pastureChange+croplandChange);
+ cropFromNatural = (croplandChange > -pastureChange) ? -(pastureChange+croplandChange) : 0;
+ pastureFromNatural = pastureChange - pastureFromCrop;
}
else {
pastureFromNatural = pastureChange;
cropFromNatural = croplandChange;
}
- // allocAllAreasForAgg(newAreaRaster, prevAreaRaster, keys, pastureFromNatural);
+ double shortfall = 0;
+ shortfall += allocAllLandCrops(newAreaRaster, keys, LandCoverType.PASTURE, LandCoverType.FOREST, pastureFromNatural * prevForestToNaturalFraction);
+ shortfall += allocAllLandCrops(newAreaRaster, keys, LandCoverType.PASTURE, LandCoverType.OTHER_NATURAL, pastureFromNatural * (1-prevForestToNaturalFraction));
+ shortfall += allocAllLandCrops(newAreaRaster, keys, LandCoverType.PASTURE, LandCoverType.CROPLAND, pastureFromCrop);
+ shortfall += allocAllLandCrops(newAreaRaster, keys, LandCoverType.CROPLAND, LandCoverType.FOREST, cropFromNatural * prevForestToNaturalFraction);
+ shortfall += allocAllLandCrops(newAreaRaster, keys, LandCoverType.CROPLAND, LandCoverType.OTHER_NATURAL, cropFromNatural * (1-prevForestToNaturalFraction));
-
- double shortfall = allocAllAreasForAgg(newAreaRaster, prevAreaRaster, keys, cropChangeTotals);
if (shortfall > 0.00001)
LogWriter.printlnError("This should never happen, due to GAMS constraint. Not able to incorporate all changes, as not enough forest or natural areas left: " + shortfall);
+
+ for (RasterKey key : keys) {
+ AreasItem newAreasRasterItem = newAreaRaster.get(key);
+ for (CropType crop : CropType.values()) {
+ newAreasRasterItem.setCropFraction(crop, newAreaAggItem.getCropFraction(crop));
+ }
+ }
}
return newAreaRaster;
}
-
- private double allocAllAreasForAgg(RasterSet<AreasItem> newAreaRaster, RasterSet<AreasItem> prevAreaRaster, Set<RasterKey> keys, CropToDouble cropChange) {
-
+
+ private double allocAllLandCrops(RasterSet<AreasItem> newAreaRaster, Set<RasterKey> keys, LandCoverType toLC, LandCoverType fromLC, double change) {
+ LandCoverType toType = toLC;
+ LandCoverType fromType = fromLC;
+
+ // reverse direction if negative
+ if (change < 0) {
+ change = -change;
+ toType = fromLC;
+ fromType = toLC;
+ }
+
double totalShortfall = 0;
Set<RasterKey> keysWithSpace = new HashSet<RasterKey>();
- CropToDouble avgCropChange = cropChange.scale(1.0/keys.size());
+ double avgChange = change/keys.size();
for (RasterKey key : keys) {
//if (DEBUG) LogWriter.println(" Processing raster key " + key);
- AreasItem prevAreasRasterItem = prevAreaRaster.get(key);
-
- AreasItem newAreasRasterItem = new AreasItem();
- double shortfall = allocRasterArea(newAreasRasterItem, prevAreasRasterItem, avgCropChange);
+ AreasItem newAreasRasterItem = newAreaRaster.get(key);
+
+ double shortfall = newAreasRasterItem.moveAreas(toType, fromType, avgChange);
if (shortfall == 0)
keysWithSpace.add(key);
else
totalShortfall += shortfall;
-
- newAreaRaster.put(key, newAreasRasterItem);
}
if (totalShortfall > 0 & keysWithSpace.size() > 0) { // more to allocate and some free areas to allocate into
- if (DEBUG) LogWriter.println("A total area of " + totalShortfall + " is unallocated, so trying again with areas with space");
- double factor = totalShortfall / cropChange.getTotal();
-
- allocAllAreasForAgg(newAreaRaster, newAreaRaster, keysWithSpace, cropChange.scale(factor));
+ if (DEBUG) LogWriter.println("A total land cover of " + totalShortfall + " is unallocated, so trying again with areas with space");
+ totalShortfall = allocAllLandCrops(newAreaRaster, keysWithSpace, toType, fromType, totalShortfall);
}
return totalShortfall;
}
- /** returns the net area of crop change that not been able to be allocated */
- private double allocRasterArea(AreasItem newAreasRasterItem, AreasItem prevAreasRasterItem, CropToDouble cropChange) {
-
- double netAddedCrop = cropChange.getTotal();
- double forestRemoved = netAddedCrop/2;
- double naturalRemoved = netAddedCrop/2;
- double prevForest = prevAreasRasterItem.getLandCoverArea(LandDataType.FOREST);
- double prevNatural = prevAreasRasterItem.getLandCoverArea(LandDataType.OTHER_NATURAL);
-
- double ratioOfChanges = 1;
-
- if (netAddedCrop > prevForest + prevNatural) {
- ratioOfChanges = (prevForest + prevNatural) / netAddedCrop;
- forestRemoved = prevForest;
- naturalRemoved = prevNatural;
- if (DEBUG) LogWriter.println("Not able to incorporate all changes for raster location: " + netAddedCrop);
- }
- else if (prevForest < forestRemoved) {
- forestRemoved = prevForest;
- naturalRemoved = netAddedCrop - forestRemoved;
- }
- else if (prevNatural < naturalRemoved) {
- naturalRemoved = prevNatural;
- forestRemoved = netAddedCrop - naturalRemoved;
- }
- else {
- // enough space to keep 50/50 allocation
- }
-
- newAreasRasterItem.setLandCoverArea(LandDataType.FOREST, prevForest - forestRemoved);
- newAreasRasterItem.setLandCoverArea(LandDataType.OTHER_NATURAL, prevNatural - naturalRemoved);
-
- for (CropType crop : CropType.values()) {
- double prevCropArea = prevAreasRasterItem.getCropArea(crop);
- Double d = cropChange.get(crop);
-
- double additionalArea = (d == null) ? 0 : d * ratioOfChanges;
- // LogWriter.println(String.format("%s, %.1f", crop, additionalArea));
-
- newAreasRasterItem.setCropArea(crop, prevCropArea + additionalArea);
- }
-
- return netAddedCrop - prevForest - prevNatural;
- }
-
private RasterSet<IntensitiesItem> allocIntensities(GamsLocationOutput gamsOutput) {
RasterSet<IntensitiesItem> newIntensityRaster = new RasterSet<IntensitiesItem>(rasterInputData.getPreviousAreas().getHeaderDetails());
@@ -201,6 +200,35 @@ public class GamsRasterOptimiser {
// as a first attempt only going to look at pasture and cereal yields, assuming other yields will be correlated to one or the other.
YieldRaster yieldRaster = rasterInputData.getYields();
RasterSet<AreasItem> cropAreaRaster = rasterInputData.getPreviousAreas();
+
+ /* checkedTotalAreas(cropAreaRaster, "before");
+ new RasterOutputer<AreasItem>(cropAreaRaster, "BeforeLandCover") {
+ @Override
+ public Double getValue(RasterKey location) {
+ AreasItem area = results.get(location);
+ if (area == null)
+ return null;
+
+ int count = 0;
+ double maxLC = 0;
+ int maxId = 0;
+ for (LandCoverType l : LandCoverType.values()) {
+ count++;
+ if (maxLC < area.getLandCoverArea(l)) {
+ maxLC = area.getLandCoverArea(l);
+ maxId = count;
+ }
+ }
+ return (double)maxId;
+ }
+
+ @Override
+ public int convertToPixelValue(double value) {
+ return (int) (1 + value * 2);
+ }
+ }.writeOutput(true); */
+
+
RasterSet<IrrigationCostItem> irrigCostRaster = rasterInputData.getIrrigationCost();
{
@@ -247,6 +275,8 @@ public class GamsRasterOptimiser {
int countFound = 0, countMissing = 0;
+
+
for (Entry<RasterKey, YieldResponsesItem> entry : yieldRaster.entrySet()) {
YieldResponsesItem yresp = entry.getValue();
RasterKey key = entry.getKey();
@@ -264,7 +294,7 @@ public class GamsRasterOptimiser {
LogWriter.printlnError("GamsRasterOptimiser: Can't find cropAreas for " + key);
continue;
}
-
+
IrrigationCostItem irrigCost = irrigCostRaster.get(key);
int cerealCat = findCategory(wheatlDivisions, yresp.getYieldMax(CropType.WHEAT));
@@ -274,7 +304,7 @@ public class GamsRasterOptimiser {
AveragingYieldResponsesItem avgYResp = aggregatedYields.lazyGet(id);
AreasItem aggAreas = aggregatedAreas.lazyGet(id);
IrrigationCostItem avgIrigCost = aggregatedIrrigCosts.lazyGet(id);
- mapping.lazyGet(id).add(entry.getKey());
+ mapping.lazyGet(id).add(key);
// Irrigation cost
double suitableAreaThisTime = cropAreas.getSuitableLand();
@@ -284,7 +314,7 @@ public class GamsRasterOptimiser {
//LogWriter.println(id + ", " + key + ", " + avgIrigCost.getIrrigCost() + ", from " + areaSoFar + ", " + suitableAreaThisTime + ", " + irrigCost.getIrrigCost());
}
- // Crop yields
+ // Crops yields and area fractions
for (CropType crop : CropType.getNonMeatTypes()) {
for (YieldType yieldType : YieldType.values()) {
avgYResp.setYield(yieldType, crop, yresp.getYield(yieldType, crop));
@@ -292,11 +322,14 @@ public class GamsRasterOptimiser {
double areaSoFar = aggAreas.getCropArea(crop);
double areaThisTime = cropAreas.getCropArea(crop);
- aggAreas.setCropArea(crop, areaSoFar + areaThisTime);
+ double updateCropland = (cropAreas.getLandCoverArea(LandCoverType.CROPLAND) + aggAreas.getLandCoverArea(LandCoverType.CROPLAND));
+
+ if (updateCropland > 0)
+ aggAreas.setCropFraction(crop, (areaSoFar + areaThisTime) / updateCropland);
}
-
- // Land cover, we only need forest and other natural
- for (LandDataType landType : Arrays.asList(LandDataType.FOREST, LandDataType.OTHER_NATURAL)) {
+
+ // Land covers ares
+ for (LandCoverType landType : LandCoverType.values()) {
double areaThisTime = cropAreas.getLandCoverArea(landType);
double areaSoFar = aggAreas.getLandCoverArea(landType);
aggAreas.setLandCoverArea(landType, areaSoFar + areaThisTime);
@@ -306,7 +339,6 @@ public class GamsRasterOptimiser {
LogWriter.println("YieldResponsesItem: " + rasterInputData.getCountryInput().getCountry() + ", countFound=" + countFound + ", countMissing=" + countMissing);
-
for (Map.Entry<Integer, Set<RasterKey>> e : mapping.entrySet()) {
LogWriter.println(e.getKey() + " category has " + e.getValue().size() + " raster areas");
if (DEBUG) for (RasterKey key : e.getValue()) {
diff --git a/src/ac/ed/lurg/landuse/AreasItem.java b/src/ac/ed/lurg/landuse/AreasItem.java
index 9b85996b..e5b6f4d7 100644
--- a/src/ac/ed/lurg/landuse/AreasItem.java
+++ b/src/ac/ed/lurg/landuse/AreasItem.java
@@ -6,7 +6,7 @@ import java.util.Map.Entry;
import ac.ed.lurg.types.CropToDouble;
import ac.ed.lurg.types.CropType;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.sac.raster.RasterItem;
/*
@@ -14,57 +14,84 @@ import ac.sac.raster.RasterItem;
*/
public class AreasItem implements RasterItem {
- private Map<CropType, Double> cropAreas = new HashMap<CropType, Double>();
- private Map<LandDataType, Double> landCoverAreas = new HashMap<LandDataType, Double>();
+ private Map<CropType, Double> cropFractions = new HashMap<CropType, Double>();
+ private Map<LandCoverType, Double> landCoverAreas = new HashMap<LandCoverType, Double>();
+
+ public void setLandCoverAreas(LandCoverItem landCover) {
+ if (landCover != null) {
+ for (LandCoverType lcType : LandCoverType.values())
+ landCoverAreas.put(lcType, landCover.getLandCover(lcType));
+ }
+ }
+
+ public void setLandCoverAreas(AreasItem otherAreasItem) {
+ for (Entry<LandCoverType, Double> entry : otherAreasItem.landCoverAreas.entrySet()) {
+ landCoverAreas.put(entry.getKey(), entry.getValue());
+ }
+ }
+
+ /** move areas from one land cover to another, return an residual not possible */
+ public double moveAreas(LandCoverType toType, LandCoverType fromType, double changeReq) {
+
+ double prevTo = getLandCoverArea(toType);
+ double prevFrom = getLandCoverArea(fromType);
+ double actuallyChanged = Math.min(prevFrom, changeReq);
+
+ setLandCoverArea(toType, prevTo + actuallyChanged);
+ setLandCoverArea(fromType, prevFrom - actuallyChanged);
+
+ return changeReq - actuallyChanged;
+ }
+
public double getCropArea(CropType c) {
- Double d = cropAreas.get(c);
+ Double d = getLandCoverArea(LandCoverType.CROPLAND) * getCropFraction(c);
+ return d;
+ }
+
+ public double getCropFraction(CropType c) {
+ Double d = cropFractions.get(c);
return d == null ? 0.0 : d;
}
+
- public void setCropArea(CropType c, double area) {
- cropAreas.put(c, area);
+ public void setCropFraction(CropType c, double area) {
+ cropFractions.put(c, area);
}
- public double getLandCoverArea(LandDataType c) {
+ public double getLandCoverArea(LandCoverType c) {
Double d = landCoverAreas.get(c);
return d == null ? 0.0 : d;
}
- public void setLandCoverArea(LandDataType c, double d) {
+ public void setLandCoverArea(LandCoverType c, double d) {
landCoverAreas.put(c, d);
}
-
- public double getTotalCropIncPastureArea() {
- double total = 0;
- for (double d : cropAreas.values()) {
- total += d;
- }
- return total;
- }
-
+
public double getSuitableLand() {
- double d = getTotalCropIncPastureArea() + getLandCoverArea(LandDataType.FOREST) + getLandCoverArea(LandDataType.OTHER_NATURAL);
+ double d = 0;
+ for (LandCoverType l : LandCoverType.values())
+ d += getLandCoverArea(l);
+
return d;
}
public double getForestToNaturalFraction() {
- double forest = getLandCoverArea(LandDataType.FOREST);
- double natural = forest + getLandCoverArea(LandDataType.OTHER_NATURAL);
+ double forest = getLandCoverArea(LandCoverType.FOREST);
+ double natural = forest + getLandCoverArea(LandCoverType.OTHER_NATURAL);
double d = forest / natural;
return d;
}
-
public CropToDouble getCropChanges(AreasItem prevAreaAggItem) {
CropToDouble changes = new CropToDouble();
- for (Entry<CropType, Double> entry : prevAreaAggItem.cropAreas.entrySet()) {
- CropType crop = entry.getKey();
- double change = getCropArea(crop) - entry.getValue();
- changes.put(crop, change);
+ for (CropType c : CropType.getCropsLessPasture()) {
+ double change = getCropArea(c) - prevAreaAggItem.getCropArea(c);
+ changes.put(c, change);
}
return changes;
}
+
}
diff --git a/src/ac/ed/lurg/landuse/LandCoverItem.java b/src/ac/ed/lurg/landuse/LandCoverItem.java
index 32eeb82c..13f6bc77 100644
--- a/src/ac/ed/lurg/landuse/LandCoverItem.java
+++ b/src/ac/ed/lurg/landuse/LandCoverItem.java
@@ -3,21 +3,21 @@ package ac.ed.lurg.landuse;
import java.util.HashMap;
import java.util.Map;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.sac.raster.RasterItem;
/** Used to hold less detailed land-cover information
* This is used in the initalisation phase, after that land-use is used */
public class LandCoverItem implements RasterItem {
- Map<LandDataType, Double> landcover = new HashMap<LandDataType, Double>();
+ Map<LandCoverType, Double> landcover = new HashMap<LandCoverType, Double>();
/** Area in Mha */
- public Double getLandCover(LandDataType landType) {
+ public Double getLandCover(LandCoverType landType) {
return landcover.get(landType);
}
- public void setLandCover(LandDataType landType, double d) {
+ public void setLandCover(LandCoverType landType, double d) {
landcover.put(landType, d);
}
diff --git a/src/ac/ed/lurg/landuse/LandCoverReader.java b/src/ac/ed/lurg/landuse/LandCoverReader.java
index 8fbd6133..11e2e761 100644
--- a/src/ac/ed/lurg/landuse/LandCoverReader.java
+++ b/src/ac/ed/lurg/landuse/LandCoverReader.java
@@ -1,14 +1,14 @@
package ac.ed.lurg.landuse;
-import ac.ed.lurg.types.LandDataType;
+import ac.ed.lurg.types.LandCoverType;
import ac.sac.raster.AbstractRasterReader;
import ac.sac.raster.RasterSet;
public class LandCoverReader extends AbstractRasterReader<LandCoverItem> {
- private LandDataType landType;
+ private LandCoverType landType;
- public LandCoverReader (RasterSet<LandCoverItem> dataset, LandDataType landType) {
+ public LandCoverReader (RasterSet<LandCoverItem> dataset, LandCoverType landType) {
super(dataset);
this.landType = landType;
}
diff --git a/src/ac/ed/lurg/types/CropToDouble.java b/src/ac/ed/lurg/types/CropToDouble.java
index b8f1176c..ada35b77 100644
--- a/src/ac/ed/lurg/types/CropToDouble.java
+++ b/src/ac/ed/lurg/types/CropToDouble.java
@@ -22,13 +22,4 @@ public class CropToDouble extends HashMap<CropType, Double> {
return total;
}
-
- public double getCroplandTotal() {
- double total = 0;
- for (CropType c : CropType.getCropsLessPasture())
- total += get(c);
-
- return total;
-
- }
}
diff --git a/src/ac/ed/lurg/types/CropType.java b/src/ac/ed/lurg/types/CropType.java
index 85c13bb3..d42471d8 100644
--- a/src/ac/ed/lurg/types/CropType.java
+++ b/src/ac/ed/lurg/types/CropType.java
@@ -40,7 +40,7 @@ public enum CropType {
Collection<CropType> comms = new HashSet<CropType>();
for (CropType c : values())
- if (c.importedCrop || !c.isMeat)
+ if (c.importedCrop && !c.isMeat)
comms.add(c);
return comms;
diff --git a/src/ac/ed/lurg/types/LandDataType.java b/src/ac/ed/lurg/types/LandCoverType.java
similarity index 79%
rename from src/ac/ed/lurg/types/LandDataType.java
rename to src/ac/ed/lurg/types/LandCoverType.java
index 9c085949..6908a9d7 100644
--- a/src/ac/ed/lurg/types/LandDataType.java
+++ b/src/ac/ed/lurg/types/LandCoverType.java
@@ -1,8 +1,7 @@
package ac.ed.lurg.types;
-public enum LandDataType {
+public enum LandCoverType {
-// LAND ("land"),
FOREST ("forest"),
OTHER_NATURAL ("otherNatural"),
CROPLAND ("cropland"),
@@ -10,7 +9,7 @@ public enum LandDataType {
private String name;
- LandDataType(String name) {
+ LandCoverType(String name) {
this.name = name;
}
diff --git a/src/ac/sac/raster/RasterSet.java b/src/ac/sac/raster/RasterSet.java
index 768c30e7..4537e751 100755
--- a/src/ac/sac/raster/RasterSet.java
+++ b/src/ac/sac/raster/RasterSet.java
@@ -5,6 +5,7 @@ import java.util.Collection;
import java.util.HashMap;
import java.util.Map;
+import ac.ed.lurg.landuse.AreasItem;
import ac.ed.lurg.utils.LogWriter;
/* Class which holds raster data. The generics used to defines the type of raster data held.
--
GitLab