no message

src/ac/ed/lurg/country/gams/GamsRasterInput.java

 package ac.ed.lurg.country.gams;
 
 import ac.ed.lurg.landuse.LandUse;
-import ac.ed.lurg.yield.YieldResponses;
+import ac.ed.lurg.yield.YieldRaster;
 import ac.sac.raster.RasterSet;
 
 public interface GamsRasterInput {
 	
-	RasterSet<YieldResponses> getYields();
+	YieldRaster getYields();
 	RasterSet<LandUse> getPreviousLandUse();
 	
 	GamsCountryInput getCountryInput();

src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java

 package ac.ed.lurg.country.gams;
 
+import java.util.ArrayList;
+import java.util.Collections;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+import java.util.Map.Entry;
+
+import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.YieldType;
+import ac.ed.lurg.utils.LogWriter;
+import ac.ed.lurg.yield.AveragingYieldResponses;
+import ac.ed.lurg.yield.YieldRaster;
+import ac.ed.lurg.yield.YieldResponses;
+import ac.sac.raster.RasterKey;
+
 public class GamsRasterOptimiser {
 	
 	private GamsRasterInput inputData;
@@ -9,9 +24,10 @@ public class GamsRasterOptimiser {
 	}
 
 	public GamsRasterOutput run() {
-		// workout similar areas
+		// workout similar areas - clustering multi-dimensional data
+		YieldRaster yieldRaster = inputData.getYields();
 		
-		// map from raster to combined location ids
+		Map<Integer, AveragingYieldResponses> mapping = findSimilarAreas(yieldRaster);
 		
 		// create GamsInput using mapping
 		
@@ -21,4 +37,71 @@ public class GamsRasterOptimiser {
 		
 		return null;
 	}
+	
+	
+	private Map<Integer, AveragingYieldResponses> findSimilarAreas(YieldRaster yieldRaster) {
+		// as a first attempt only going to look at pasture and cereal yields, assuming other yields will be correlated to one or the other.
+		List<Double> cerealDivisions = getDivisions(yieldRaster, CropType.CEREALS, 2);
+		List<Double> pastureDivisions = getDivisions(yieldRaster, CropType.MEAT_OR_PASTURE, 2);
+		
+		int numCerealCats = (1+cerealDivisions.size());
+		LogWriter.println("Found " + numCerealCats * (1+pastureDivisions.size()) + " categories");
+		
+		Map<Integer, AveragingYieldResponses> aggregatedData = new HashMap<Integer, AveragingYieldResponses>();
+		
+		for (Entry<RasterKey, YieldResponses> entry : yieldRaster.entrySet()) {
+			YieldResponses yresp = entry.getValue();
+		
+			int cerealCat = findCategory(cerealDivisions, yresp.getYieldMax(CropType.CEREALS));
+			int pastureCat = findCategory(pastureDivisions, yresp.getYieldMax(CropType.MEAT_OR_PASTURE));
+			int id = cerealCat + pastureCat * numCerealCats;
+
+			AveragingYieldResponses avgYResp = aggregatedData.get(id);
+			if (avgYResp == null) {
+				avgYResp = new AveragingYieldResponses();
+				aggregatedData.put(id, avgYResp);	
+			}
+			
+			avgYResp.addMapping(entry.getKey());
+
+			for (CropType crop : CropType.values()) {
+				for (YieldType yieldType : YieldType.values()) {
+					avgYResp.setYield(yieldType, crop, yresp.getYield(yieldType, crop));
+				}
+			}
+		}
+		
+		return aggregatedData;
+	}
+
+	private int findCategory(List<Double> divisions, double yield) {
+		int category;
+		int numDivisions = divisions.size();
+		
+		for (category = 0; category<numDivisions; category++) {
+			if (yield > divisions.get(category)) {
+				break;
+			}
+		}
+		return category;
+	}
+	
+	private List<Double> getDivisions(YieldRaster yieldRaster, CropType crop, int numCategories) {
+		List<Double> yieldValues = new ArrayList<Double>();
+
+		for (YieldResponses yresp : yieldRaster.values()) {
+			yieldValues.add(yresp.getYieldMax(crop));
+		}
+		
+		Collections.sort(yieldValues);
+
+		List<Double> divisions = new ArrayList<Double>();
+		
+		for (int i=1; i < numCategories; i++) {
+			double d = yieldValues.get(yieldValues.size()/numCategories);
+			divisions.add(d);
+		}	
+		
+		return divisions;
+	}
 }

src/ac/ed/lurg/yield/AveragingYieldResponse.java

+package ac.ed.lurg.yield;
+
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.types.YieldType;
+
+public class AveragingYieldResponse extends YieldResponse {
+	private Map<YieldType, Integer> numAreas = new HashMap<YieldType, Integer>();
+
+	// we are assuming that the areas are all the same size, i.e. countries boundaries follow grid structure
+	public void setYield(YieldType type, double yield) {
+		double previousYield = yields.get(type);
+		int prevNumAreas = numAreas.get(type);
+
+		double newAvgYield = (previousYield*prevNumAreas + yield)/(prevNumAreas+1);
+		
+		yields.put(type, newAvgYield);
+		numAreas.put(type, prevNumAreas+1);
+	}
+}

src/ac/ed/lurg/yield/AveragingYieldResponses.java

+package ac.ed.lurg.yield;
+
+import java.util.Collection;
+import java.util.HashSet;
+
+import ac.ed.lurg.types.CropType;
+import ac.sac.raster.RasterKey;
+
+public class AveragingYieldResponses extends YieldResponses {
+
+	private Collection<RasterKey> rasterKeys = new HashSet<RasterKey>();
+	
+	@Override
+	YieldResponse getYieldResponseForCrop(CropType crop) {
+		YieldResponse yresp = yieldResponses.get(crop);
+		if (yresp == null) {
+			yresp = new AveragingYieldResponse();
+			yieldResponses.put(crop, yresp);
+		}
+		return yresp;
+	}
+
+	public void addMapping(RasterKey key) {
+		rasterKeys.add(key);
+	}
+	
+	public Collection<RasterKey> getRasterKeys() {
+		return rasterKeys;
+	}
+}

src/ac/ed/lurg/yield/YieldRaster.java

+package ac.ed.lurg.yield;
+
+import ac.ed.lurg.country.Country;
+import ac.sac.raster.RasterHeaderDetails;
+import ac.sac.raster.RasterSet;
+
+public class YieldRaster extends RasterSet<YieldResponses> {
+
+	private static final long serialVersionUID = -4483319360580776344L;
+	
+	public YieldRaster(RasterHeaderDetails header) {
+		super(header);
+	}
+	
+	protected YieldResponses createRasterData() {
+		return new YieldResponses();
+	}
+	
+	public YieldRaster getRasterForCountry(Country c) {
+		// TODO
+		return null;
+	}
+}

src/ac/ed/lurg/yield/YieldResponse.java

@@ -6,7 +6,7 @@ import java.util.Map;
 import ac.ed.lurg.types.YieldType;
 
 public class YieldResponse {
-	private Map<YieldType, Double> yields = new HashMap<YieldType, Double>();
+	Map<YieldType, Double> yields = new HashMap<YieldType, Double>();
 
 	public void setYield(YieldType type, double yield) {
 		yields.put(type, yield);

src/ac/ed/lurg/yield/YieldResponseMapReader.java

@@ -22,13 +22,7 @@ public class YieldResponseMapReader extends AbstractRasterReader<YieldResponses>
 		if (DEBUG) LogWriter.println("Creating RasterDataset col:"+header.getNcolumns() + ", rows:" + header.getNrows());
 
 		if (dataset == null) {
-			dataset = new RasterSet<YieldResponses>(header) {
-				private static final long serialVersionUID = -202195565674146064L;
-
-				protected YieldResponses createRasterData() {
-					return new YieldResponses();
-				}
-			};
+			dataset = new YieldRaster(header);
 		}
 	}
 

src/ac/ed/lurg/yield/YieldResponses.java

@@ -8,13 +8,13 @@ import ac.ed.lurg.types.YieldType;
 import ac.sac.raster.RasterItem;
 
 public class YieldResponses implements RasterItem {
-	private Map<CropType, YieldResponse> yieldResponses = new HashMap<CropType, YieldResponse>();
+	Map<CropType, YieldResponse> yieldResponses = new HashMap<CropType, YieldResponse>();
 
 	public void setYield(YieldType type, CropType crop, double yield) {
 		getYieldResponseForCrop(crop).setYield(type, yield);
 	}
 
-	private YieldResponse getYieldResponseForCrop(CropType crop) {
+	YieldResponse getYieldResponseForCrop(CropType crop) {
 		YieldResponse yresp = yieldResponses.get(crop);
 		if (yresp == null) {
 			yresp = new YieldResponse();
@@ -23,20 +23,24 @@ public class YieldResponses implements RasterItem {
 		return yresp;
 	}
 	
+	public double getYield(YieldType yieldType, CropType crop) {
+		return getYieldResponseForCrop(crop).getYield(yieldType);
+	}	
+
 	public double getYieldNone(CropType crop) {
-		return getYieldResponseForCrop(crop).getYield(YieldType.NONE);
+		return getYield(YieldType.NONE, crop);
 	}	
 	
 	public double getYieldFertOnly(CropType crop) {
-		return getYieldResponseForCrop(crop).getYield(YieldType.FERT_MAX_ONLY);
+		return getYield(YieldType.FERT_MAX_ONLY, crop);
 	}
 	
 	public double getYieldIrrigOnly(CropType crop) {
-		return getYieldResponseForCrop(crop).getYield(YieldType.IRRIG_MAX_ONLY);
+		return getYield(YieldType.IRRIG_MAX_ONLY, crop);
 	}
 	
 	public double getYieldMax(CropType crop) {
-		return getYieldResponseForCrop(crop).getYield(YieldType.FERT_IRRIG_MAX);
+		return getYield(YieldType.FERT_IRRIG_MAX, crop);
 	}
 	
 	public double getFertParam(CropType crop) {