diff --git a/src/ac/ed/lurg/country/gams/GamsRasterInput.java b/src/ac/ed/lurg/country/gams/GamsRasterInput.java
index 8ca8a301bf7c9f9a3578d02c02cd130215f837fb..1292d3b0a654dedcd49baaa517042939e0f5deb9 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterInput.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterInput.java
@@ -1,12 +1,12 @@
 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();
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
index 0d573b83b194d2b1411b80d012c88b80bc9cd7e6..74877b74c9b895c33e024035c35f8dc50f36ae24 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
@@ -1,5 +1,20 @@
 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;
+	}
 }
diff --git a/src/ac/ed/lurg/yield/AveragingYieldResponse.java b/src/ac/ed/lurg/yield/AveragingYieldResponse.java
new file mode 100644
index 0000000000000000000000000000000000000000..7f35276bc4792b85dc1bd1b935fd0517fea3dd6d
--- /dev/null
+++ b/src/ac/ed/lurg/yield/AveragingYieldResponse.java
@@ -0,0 +1,21 @@
+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);
+	}
+}
diff --git a/src/ac/ed/lurg/yield/AveragingYieldResponses.java b/src/ac/ed/lurg/yield/AveragingYieldResponses.java
new file mode 100644
index 0000000000000000000000000000000000000000..9d34d694db8c3ca1191aa7ec6679fe924533e1c0
--- /dev/null
+++ b/src/ac/ed/lurg/yield/AveragingYieldResponses.java
@@ -0,0 +1,30 @@
+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;
+	}
+}
diff --git a/src/ac/ed/lurg/yield/YieldRaster.java b/src/ac/ed/lurg/yield/YieldRaster.java
new file mode 100644
index 0000000000000000000000000000000000000000..2642384190469d19500ccb5a49c784847464b621
--- /dev/null
+++ b/src/ac/ed/lurg/yield/YieldRaster.java
@@ -0,0 +1,23 @@
+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;
+	}
+}
diff --git a/src/ac/ed/lurg/yield/YieldResponse.java b/src/ac/ed/lurg/yield/YieldResponse.java
index 6f621e093ff3e32e6233f5f7f57b4a39aa43bd38..bec1c20f670613e0c4ad5d008d947f48a681d74a 100644
--- a/src/ac/ed/lurg/yield/YieldResponse.java
+++ b/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);
diff --git a/src/ac/ed/lurg/yield/YieldResponseMapReader.java b/src/ac/ed/lurg/yield/YieldResponseMapReader.java
index 932f3758b8fa6c4798d95b817ad4ba4fc149c542..6a19695eb886b98c197b7a38225a8db4b7caa965 100644
--- a/src/ac/ed/lurg/yield/YieldResponseMapReader.java
+++ b/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);
 		}
 	}
 
diff --git a/src/ac/ed/lurg/yield/YieldResponses.java b/src/ac/ed/lurg/yield/YieldResponses.java
index 9f19a38cbb321d43353c38bc0ccad6c4328b2b84..799ec7a6cd71581652e5deba0fc533315af5ca85 100644
--- a/src/ac/ed/lurg/yield/YieldResponses.java
+++ b/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) {