From fc35df6b1429333e6dff3c348b898e124c3bec9e Mon Sep 17 00:00:00 2001
From: Peter Alexander <peter@blackhillock.co.uk>
Date: Tue, 9 Dec 2014 19:25:39 +0000
Subject: [PATCH] no message

---
 GAMS/IntExtOpt.gms                            |  36 ++---
 src/ac/ed/lurg/ModelConfig.java               |   1 +
 src/ac/ed/lurg/ModelMain.java                 |  21 ++-
 .../ed/lurg/country/gams/GamsInputData.java   |   7 +-
 .../ed/lurg/country/gams/GamsOutputData.java  |  15 ++
 src/ac/ed/lurg/landuse/LandUseDataPoint.java  |  72 ++++++++++
 src/ac/ed/lurg/landuse/LandUseItem.java       |  50 +++++++
 src/ac/ed/lurg/landuse/LandUseSet.java        |  12 ++
 src/ac/ed/lurg/types/YieldType.java           |  21 +++
 src/ac/ed/lurg/yield/YieldResponse.java       |  37 +++++
 .../ed/lurg/yield/YieldResponseMapReader.java |  37 +++++
 src/ac/sac/raster/AbstractRasterReader.java   | 136 ++++++++++++++++++
 src/ac/sac/raster/RasterAreaCircle.java       |  17 +++
 src/ac/sac/raster/RasterHeaderDetails.java    | 117 +++++++++++++++
 src/ac/sac/raster/RasterItem.java             |   7 +
 src/ac/sac/raster/RasterKey.java              |  69 +++++++++
 src/ac/sac/raster/RasterSet.java              | 113 +++++++++++++++
 17 files changed, 747 insertions(+), 21 deletions(-)
 create mode 100644 src/ac/ed/lurg/country/gams/GamsOutputData.java
 create mode 100644 src/ac/ed/lurg/landuse/LandUseDataPoint.java
 create mode 100644 src/ac/ed/lurg/landuse/LandUseItem.java
 create mode 100644 src/ac/ed/lurg/landuse/LandUseSet.java
 create mode 100644 src/ac/ed/lurg/types/YieldType.java
 create mode 100644 src/ac/ed/lurg/yield/YieldResponse.java
 create mode 100644 src/ac/ed/lurg/yield/YieldResponseMapReader.java
 create mode 100755 src/ac/sac/raster/AbstractRasterReader.java
 create mode 100755 src/ac/sac/raster/RasterAreaCircle.java
 create mode 100755 src/ac/sac/raster/RasterHeaderDetails.java
 create mode 100755 src/ac/sac/raster/RasterItem.java
 create mode 100755 src/ac/sac/raster/RasterKey.java
 create mode 100755 src/ac/sac/raster/RasterSet.java

diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index a6c3795e..35bdcdc5 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -44,30 +44,30 @@ $gdxin
        fertI(crop, location)              fertilizer intensity for each crop - factor between 0 and 1
        irrigI(crop, location)             irrigation intensity for each crop - factor between 0 and 1
        otherIntensity(crop, location)     other intensity for each crop - unit less
-       feedAmount(crop)         amount of feed for each crop - t
-       netImportAmount(crop)    export of crop - t
+       feedAmount(crop)                   amount of feed for each crop - t
+       netImportAmount(crop)              export of crop - t
        yield(crop, location)              yield per area for each crop - t per ha
        unit_energy(crop, location)        energy per area for each crop - energy
-       energy                   total input energy - energy;
+       energy                             total input energy - energy;
  
  POSITIVE VARIABLE area, fertI, irrigI, otherIntensity, feedAmount;
  
  EQUATIONS
-       UNIT_ENERGY_EQ(crop, location)                       energy per area
-       YIELD_EQ(crop, location)                             yield given chosen intensity
-       CROP_DEMAND_CONSTRAINT(crop_less_pasture)  satisfy demand for crop 
-       MEAT_DEMAND_CONSTRAINT                     satisfy demand for crop
-       MAX_FERT_INTENSITY_CONSTRAINT(crop, location)        constraint on maximum fertilizer intensity
-       MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location)       constraint on maximum irrigation intensity
-       CROPLAND_CHANGE_CONSTRAINT(location)                 constraint on the rate of land use change 
-       PASTURE_CHANGE_CONSTRAINT(location)                  constraint on the rate of land use change 
-       AGRI_LAND_INCREASE_CONSTRAINT(location)                constraint expansion of agricultural land
-       AGRI_LAND_DECREASE_CONSTRAINT(location)                constraint expansion of agricultural land
-       NON_FEED_CROP_CONSTRAINT(not_feed_crops)   constraint to set non feed crop feed usage to zero
-       MAX_NET_IMPORT_CONSTRAINT(crop)            constraint on max net imports
-       MIN_NET_IMPORT_CONSTRAINT(crop)            constraint on min net imports
-       MIN_FEED_CONSTRAINT                        constraint on min feed rate
-       ENERGY_EQ                                  total energy objective function;
+       UNIT_ENERGY_EQ(crop, location)                   energy per area
+       YIELD_EQ(crop, location)                         yield given chosen intensity
+       CROP_DEMAND_CONSTRAINT(crop_less_pasture)        satisfy demand for crop 
+       MEAT_DEMAND_CONSTRAINT                           satisfy demand for crop
+       MAX_FERT_INTENSITY_CONSTRAINT(crop, location)    constraint on maximum fertilizer intensity
+       MAX_IRRIG_INTENSITY_CONSTRAINT(crop, location)   constraint on maximum irrigation intensity
+       CROPLAND_CHANGE_CONSTRAINT(location)             constraint on the rate of land use change 
+       PASTURE_CHANGE_CONSTRAINT(location)              constraint on the rate of land use change 
+       AGRI_LAND_INCREASE_CONSTRAINT(location)          constraint expansion of agricultural land
+       AGRI_LAND_DECREASE_CONSTRAINT(location)          constraint expansion of agricultural land
+       NON_FEED_CROP_CONSTRAINT(not_feed_crops)         constraint to set non feed crop feed usage to zero
+       MAX_NET_IMPORT_CONSTRAINT(crop)                  constraint on max net imports
+       MIN_NET_IMPORT_CONSTRAINT(crop)                  constraint on min net imports
+       MIN_FEED_CONSTRAINT                              constraint on min feed rate
+       ENERGY_EQ                                        total energy objective function;
  
  UNIT_ENERGY_EQ(crop, location) .. unit_energy(crop, location) =E= power((1+ fertI(crop, location)+irrigI(crop, location)+otherIntensity(crop, location)), 2);
  
diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java
index 1208fb42..1894eb05 100644
--- a/src/ac/ed/lurg/ModelConfig.java
+++ b/src/ac/ed/lurg/ModelConfig.java
@@ -90,6 +90,7 @@ public class ModelConfig {
 	public static final String DEMAND_CURVES_FILE = DATA_DIR + File.separator + "com_curves.csv";
 	public static final String SSP_FILE = DATA_DIR + File.separator + "ssp.csv";
 	public static final String COUNTRY_DATA_FILE = DATA_DIR + File.separator + "country_data.csv";
+	public static final String YIELD_DIR = DATA_DIR + File.separator + "yields";
 	public static final String REF_YIELD_FILE = DATA_DIR + File.separator + "ref_yields.csv";
 	public static final String BASELINE_CONSUMP_FILE = DATA_DIR + File.separator + "base_consump.csv";
 
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 621f1558..c9e9c175 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -1,12 +1,17 @@
 package ac.ed.lurg;
 
+import java.io.File;
 import java.util.Collection;
 
 import ac.ed.lurg.country.CountryAgent;
 import ac.ed.lurg.country.CountryAgentCreator;
 import ac.ed.lurg.demand.DemandManager;
 import ac.ed.lurg.types.ModelFitType;
+import ac.ed.lurg.types.YieldType;
 import ac.ed.lurg.yield.YieldManager;
+import ac.ed.lurg.yield.YieldResponse;
+import ac.ed.lurg.yield.YieldResponseMapReader;
+import ac.sac.raster.RasterSet;
 
 public class ModelMain {
 
@@ -37,10 +42,24 @@ public class ModelMain {
 	}
 
 	private void doTimestep(int timestep) {
+		RasterSet<YieldResponse> yieldSurfaces = getYieldSurfaces(timestep);
+		
 		for (CountryAgent ca : countryAgents) {
 			ca.determineProduction(timestep);
 		}
 		
 		// examine global trade balance
-	}	
+	}
+
+	private RasterSet<YieldResponse> getYieldSurfaces(int timestep) {
+		String rootDir = ModelConfig.YIELD_DIR + File.separator + (timestep + ModelConfig.START_TIMESTEP) + File.separator;
+		RasterSet<YieldResponse> yieldSurfaces = null;
+		
+		for (YieldType yieldType : YieldType.values()) {
+			YieldResponseMapReader yieldReader = new YieldResponseMapReader(yieldSurfaces, yieldType);
+			yieldSurfaces = yieldReader.getRasterDataFromFile(rootDir + yieldType.getFileName()); 
+		}
+		
+		return yieldSurfaces;
+	}
 }
diff --git a/src/ac/ed/lurg/country/gams/GamsInputData.java b/src/ac/ed/lurg/country/gams/GamsInputData.java
index 4ea87e1a..2de64a49 100644
--- a/src/ac/ed/lurg/country/gams/GamsInputData.java
+++ b/src/ac/ed/lurg/country/gams/GamsInputData.java
@@ -10,16 +10,19 @@ public interface GamsInputData {
 	Map<Integer, Map<CropType, Double>> getYieldFertOnly();
 	Map<Integer, Map<CropType, Double>> getYieldIrrigOnly();
 	Map<Integer, Map<CropType, Double>> getYieldBoth();
-	Map<Integer, Map<CropType, Double>> getPreviousCropArea();
-	
 	Map<Integer, Map<CropType, Double>> getFertParam();
 	Map<Integer, Map<CropType, Double>> getIrrigParam();
+
+	Map<Integer, Map<CropType, Double>> getPreviousCropArea();
+	
 	
 	Map<CropType, Double> getProjectedDemand();
 	Map<CropType, Double> getWorldInputEnergy();
 	Map<CropType, Double> getMaxNetImport();
 	Map<CropType, Double> getMinNetImport();
 	
+	// limits to areas for each location and crop
+	
 	double getMeatEfficiency();
 	double getMaxLandUseChange();
 	double getMaxIntensity();
diff --git a/src/ac/ed/lurg/country/gams/GamsOutputData.java b/src/ac/ed/lurg/country/gams/GamsOutputData.java
new file mode 100644
index 00000000..dc58d5f5
--- /dev/null
+++ b/src/ac/ed/lurg/country/gams/GamsOutputData.java
@@ -0,0 +1,15 @@
+package ac.ed.lurg.country.gams;
+
+import java.util.Map;
+
+import ac.ed.lurg.landuse.LandUseItem;
+import ac.ed.lurg.types.CropType;
+
+public class GamsOutputData {
+	int status;
+	
+	Map<Integer, LandUseItem> landuses;  // data mapped from id (not raster)
+	
+	Map<CropType, Double> feedAmounts;
+	Map<CropType, Double> netImports;
+}
diff --git a/src/ac/ed/lurg/landuse/LandUseDataPoint.java b/src/ac/ed/lurg/landuse/LandUseDataPoint.java
new file mode 100644
index 00000000..8e093bd5
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/LandUseDataPoint.java
@@ -0,0 +1,72 @@
+package ac.ed.lurg.landuse;
+
+public class LandUseDataPoint {
+	double area;
+	double fertiliserIntensity;
+	double irrigationIntensity;
+	double otherIntensity;
+	
+	public LandUseDataPoint(double area, double fertiliserIntensity, double irrigationIntensity, double otherIntensity) {
+		super();
+		this.area = area;
+		this.fertiliserIntensity = fertiliserIntensity;
+		this.irrigationIntensity = irrigationIntensity;
+		this.otherIntensity = otherIntensity;
+	}
+
+	public double getArea() {
+		return area;
+	}
+
+	public double getFertiliserIntensity() {
+		return fertiliserIntensity;
+	}
+
+	public double getIrrigationIntensity() {
+		return irrigationIntensity;
+	}
+
+	public double getOtherIntensity() {
+		return otherIntensity;
+	}
+
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		long temp;
+		temp = Double.doubleToLongBits(area);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		temp = Double.doubleToLongBits(fertiliserIntensity);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		temp = Double.doubleToLongBits(irrigationIntensity);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		temp = Double.doubleToLongBits(otherIntensity);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		return result;
+	}
+
+	@Override
+	public boolean equals(Object obj) {
+		if (this == obj)
+			return true;
+		if (obj == null)
+			return false;
+		if (getClass() != obj.getClass())
+			return false;
+		LandUseDataPoint other = (LandUseDataPoint) obj;
+		if (Double.doubleToLongBits(area) != Double
+				.doubleToLongBits(other.area))
+			return false;
+		if (Double.doubleToLongBits(fertiliserIntensity) != Double
+				.doubleToLongBits(other.fertiliserIntensity))
+			return false;
+		if (Double.doubleToLongBits(irrigationIntensity) != Double
+				.doubleToLongBits(other.irrigationIntensity))
+			return false;
+		if (Double.doubleToLongBits(otherIntensity) != Double
+				.doubleToLongBits(other.otherIntensity))
+			return false;
+		return true;
+	}
+}
diff --git a/src/ac/ed/lurg/landuse/LandUseItem.java b/src/ac/ed/lurg/landuse/LandUseItem.java
new file mode 100644
index 00000000..51db7b62
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/LandUseItem.java
@@ -0,0 +1,50 @@
+package ac.ed.lurg.landuse;
+
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.types.CropType;
+import ac.sac.raster.RasterItem;
+
+public class LandUseItem implements RasterItem {
+
+	Map<CropType, LandUseDataPoint> landUses = new HashMap<CropType, LandUseDataPoint>();
+	
+	public LandUseDataPoint getLandUses(CropType crop) {
+		return landUses.get(crop);
+	}
+	
+	public void setLandUses(CropType crop, LandUseDataPoint landUse) {
+		landUses.put(crop, landUse);
+	}
+
+	public void setLandUses(Map<CropType, LandUseDataPoint> landUses) {
+		this.landUses = landUses;
+	}
+
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		result = prime * result
+				+ ((landUses == null) ? 0 : landUses.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;
+		LandUseItem other = (LandUseItem) obj;
+		if (landUses == null) {
+			if (other.landUses != null)
+				return false;
+		} else if (!landUses.equals(other.landUses))
+			return false;
+		return true;
+	}
+}
diff --git a/src/ac/ed/lurg/landuse/LandUseSet.java b/src/ac/ed/lurg/landuse/LandUseSet.java
new file mode 100644
index 00000000..192d0b52
--- /dev/null
+++ b/src/ac/ed/lurg/landuse/LandUseSet.java
@@ -0,0 +1,12 @@
+package ac.ed.lurg.landuse;
+
+import ac.sac.raster.RasterSet;
+
+public class LandUseSet extends RasterSet<LandUseItem> {
+
+	private static final long serialVersionUID = -2682317143041974233L;
+
+	protected LandUseItem createRasterData() {
+		return new LandUseItem();
+	}
+}
diff --git a/src/ac/ed/lurg/types/YieldType.java b/src/ac/ed/lurg/types/YieldType.java
new file mode 100644
index 00000000..4333084a
--- /dev/null
+++ b/src/ac/ed/lurg/types/YieldType.java
@@ -0,0 +1,21 @@
+package ac.ed.lurg.types;
+
+
+public enum YieldType {
+	NONE("y_0_0.asc"),
+	FERT_MID_ONLY("y_mid_0.asc"),
+	FERT_MAX_ONLY("y_max_0.asc"),
+	IRRIG_MID_ONLY("y_0_mid.asc"),
+	IRRIG_MAX_ONLY("y_0_max.asc"),
+	FERT_IRRIG_MAX("y_max_max.asc");
+	
+	private String fileName;
+
+	YieldType(String fileName) {
+		this.fileName = fileName;
+	}
+
+	public String getFileName() {
+		return fileName;
+	}
+}
diff --git a/src/ac/ed/lurg/yield/YieldResponse.java b/src/ac/ed/lurg/yield/YieldResponse.java
new file mode 100644
index 00000000..88fcd9bb
--- /dev/null
+++ b/src/ac/ed/lurg/yield/YieldResponse.java
@@ -0,0 +1,37 @@
+package ac.ed.lurg.yield;
+
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.types.YieldType;
+import ac.sac.raster.RasterItem;
+
+public class YieldResponse implements RasterItem {
+	private Map<YieldType, Double> yields = new HashMap<YieldType, Double>();
+
+	public void setYield(YieldType type, double yield) {
+		yields.put(type, yield);
+	}
+	
+	public double getYieldNone() {
+		return yields.get(YieldType.NONE);
+	}	
+	public double getYieldFertOnly() {
+		return yields.get(YieldType.FERT_MAX_ONLY);
+	}
+	public double getYieldIrrigOnly() {
+		return yields.get(YieldType.IRRIG_MAX_ONLY);
+	}
+	public double getYieldMax() {
+		return yields.get(YieldType.FERT_IRRIG_MAX);
+	}
+	
+	public double getFertParam() {
+		double fertMin = 0, fertMid = 0.5, fertMax = 1;
+		return Math.log((yields.get(YieldType.FERT_MID_ONLY)-getYieldNone())/(getYieldFertOnly()-getYieldNone()))/Math.log((fertMid - fertMin)/(fertMax - fertMin))/fertMax;
+	}
+	public double getIrrigParam() {
+		double irrigMin = 0, irrigMid = 0.5, irrigMax = 1;
+		return Math.log((yields.get(YieldType.IRRIG_MID_ONLY)-getYieldNone())/(getYieldIrrigOnly()-getYieldNone()))/Math.log((irrigMid - irrigMin)/(irrigMax - irrigMin))/irrigMax;
+	}
+}
diff --git a/src/ac/ed/lurg/yield/YieldResponseMapReader.java b/src/ac/ed/lurg/yield/YieldResponseMapReader.java
new file mode 100644
index 00000000..2fa41a4d
--- /dev/null
+++ b/src/ac/ed/lurg/yield/YieldResponseMapReader.java
@@ -0,0 +1,37 @@
+package ac.ed.lurg.yield;
+
+import ac.ed.lurg.types.YieldType;
+import ac.ed.lurg.utils.LogWriter;
+import ac.sac.raster.AbstractRasterReader;
+import ac.sac.raster.RasterHeaderDetails;
+import ac.sac.raster.RasterSet;
+
+public class YieldResponseMapReader extends AbstractRasterReader<YieldResponse> {	
+
+	private YieldType yieldType;
+	
+	public YieldResponseMapReader (RasterSet<YieldResponse> dataset, YieldType yieldType) {
+		super(dataset);
+		this.yieldType = yieldType;
+	}
+
+	protected void createDataSet(RasterHeaderDetails header) {
+		if (DEBUG) LogWriter.println("Creating RasterDataset col:"+header.getNcolumns() + ", rows:" + header.getNrows());
+
+		if (dataset == null) {
+			dataset = new RasterSet<YieldResponse>(header) {
+				private static final long serialVersionUID = -202195565674146064L;
+
+				protected YieldResponse createRasterData() {
+					return new YieldResponse();
+				}
+			};
+		}
+	}
+
+	@Override
+	public void setData(YieldResponse item, String token) {
+		double d = Double.parseDouble(token);
+		item.setYield(yieldType, d);
+	}
+}
\ No newline at end of file
diff --git a/src/ac/sac/raster/AbstractRasterReader.java b/src/ac/sac/raster/AbstractRasterReader.java
new file mode 100755
index 00000000..284821a6
--- /dev/null
+++ b/src/ac/sac/raster/AbstractRasterReader.java
@@ -0,0 +1,136 @@
+package ac.sac.raster;
+
+import java.io.BufferedReader;
+import java.io.FileReader;
+import java.io.IOException;
+import java.util.Collection;
+
+import ac.ed.lurg.utils.LogWriter;
+
+public abstract class AbstractRasterReader<D extends RasterItem> {
+	protected static final boolean DEBUG = false;
+	private static final int HEADER_LENGTH = 6;
+	private static final int HEADER_DATA_COL = 1;
+	private static final int NCOLUMN_ROW = 0;
+	private static final int NROW_ROW = 1;
+	private static final int XLL_CORNER = 2;
+	private static final int YLL_CORNER = 3;
+	private static final int CELL_SIZE = 4;
+	private static final int NODATA_VALUE_ROW = 5;
+			
+	protected RasterSet<D> dataset;
+	
+	public AbstractRasterReader () {
+	}
+
+	public AbstractRasterReader (RasterSet<D> dataset) {
+		this.dataset = dataset;
+	}
+	
+	protected String[] parseLine(String line) {
+		return line.split(" +|,");
+	}
+	
+	protected RasterHeaderDetails handleHeader(BufferedReader br) throws IOException {
+		String headerLine = new String();
+		int ncolumns=0;
+		int nrows=0;
+		double xllCorner=0.0;
+		double yllCorner = 0.0;
+		int cellSize = 1000;
+		String nodataString=null;
+		
+		for (int i=0; i<HEADER_LENGTH; i++) {
+			headerLine = br.readLine();
+			if (headerLine == null)
+				throw new IOException("Raster file format problem. Too short, not even a full header:" + i);
+			
+			String[] tokens = parseLine(headerLine);
+			
+			if (DEBUG)
+				for (String token : tokens)
+					LogWriter.println(token);		
+
+			switch (i) {
+			case NCOLUMN_ROW:
+				ncolumns = Integer.parseInt(tokens[HEADER_DATA_COL]);
+				break;
+			case NROW_ROW:
+				nrows = Integer.parseInt(tokens[HEADER_DATA_COL]);
+				break;
+			case XLL_CORNER:
+				xllCorner = Double.parseDouble(tokens[HEADER_DATA_COL]);
+				break;
+			case YLL_CORNER:
+				yllCorner = Double.parseDouble(tokens[HEADER_DATA_COL]);
+				break;
+			case CELL_SIZE:
+				cellSize =  Integer.parseInt(tokens[HEADER_DATA_COL]);
+				break;
+			case NODATA_VALUE_ROW:
+				nodataString = tokens[HEADER_DATA_COL];
+				break;
+			}
+		}
+		
+		if (DEBUG) LogWriter.println("Creating RasterDataset col:"+ncolumns + ", rows:" + nrows);
+		
+		RasterHeaderDetails headerDetails = new RasterHeaderDetails(ncolumns, nrows, xllCorner, yllCorner, cellSize, nodataString);
+		createDataSet(headerDetails);
+		return headerDetails;
+	}
+	
+	protected abstract void createDataSet(RasterHeaderDetails headerDetails);
+	
+	public RasterSet<D> getRasterDataFromFile(String filename) {
+
+		long startTime = System.currentTimeMillis();
+		int row = 0, col = 0;
+		
+		try {
+			BufferedReader in = new BufferedReader(new FileReader(filename));
+			RasterHeaderDetails header = handleHeader(in);
+					
+			String line;
+			
+			while ((line=in.readLine()) != null) {
+				String[] tokens = parseLine(line);
+
+				for (String token : tokens) {
+					
+					if (dataset.getHeaderDetails().getNodataString().equals(token)) {
+						if (DEBUG) LogWriter.println(",");
+					}
+					else {
+						if (DEBUG) LogWriter.println(token + ",");
+
+						Collection <D> dataPoints = dataset.getCollection(col, row, header);
+						for (D data : dataPoints) {
+							setData(data, token);
+						}
+					}					
+					col++;
+				}
+				row++;
+				if (DEBUG) LogWriter.println(filename + ", row:" + row + ", col:" + col);
+				col=0;
+			}
+			in.close();
+		}
+		catch (Exception e) {
+			LogWriter.println("row:" + row + ", col:" + col);
+			e.printStackTrace();
+		}
+
+		LogWriter.println("Reading " + filename + ", took " + (System.currentTimeMillis() - startTime) + " ms");
+
+		return dataset;
+	}
+	
+	protected boolean areHeadersConsistent(RasterHeaderDetails header) {
+		return dataset.isConstistentWithHeader(header);
+	}
+	
+	public abstract void setData(D data, String value);
+
+}
diff --git a/src/ac/sac/raster/RasterAreaCircle.java b/src/ac/sac/raster/RasterAreaCircle.java
new file mode 100755
index 00000000..3d38892f
--- /dev/null
+++ b/src/ac/sac/raster/RasterAreaCircle.java
@@ -0,0 +1,17 @@
+package ac.sac.raster;
+
+public class RasterAreaCircle {
+
+	private double radius;
+	private RasterKey center;
+	
+	public RasterAreaCircle (double radius, RasterKey center) {
+		this.radius = radius;
+		this.center = center;
+	}
+
+	public boolean contains(RasterKey location) {
+		double distance = center.getDistanceTo(location);
+		return distance<=radius;
+	}
+}
diff --git a/src/ac/sac/raster/RasterHeaderDetails.java b/src/ac/sac/raster/RasterHeaderDetails.java
new file mode 100755
index 00000000..c45b0a52
--- /dev/null
+++ b/src/ac/sac/raster/RasterHeaderDetails.java
@@ -0,0 +1,117 @@
+package ac.sac.raster;
+
+public class RasterHeaderDetails {
+	private int ncolumns;
+	private int nrows;
+	private double xllCorner;
+	private double yllCorner;
+	private int cellSize;
+	private String nodataString;
+	
+	public RasterHeaderDetails(int ncolumns, int nrows, double xllCorner, double yllCorner, int cellSize, String nodataString) {
+		super();
+		this.ncolumns = ncolumns;
+		this.nrows = nrows;
+		this.xllCorner = xllCorner;
+		this.yllCorner = yllCorner;
+		this.cellSize = cellSize;
+		this.nodataString = nodataString;
+	}
+
+	public RasterHeaderDetails() {
+	}
+	
+	public RasterHeaderDetails copy() {
+		return new RasterHeaderDetails(ncolumns, nrows, xllCorner, yllCorner, cellSize, nodataString);
+	}
+
+	public void incrementNcolumns() {
+		ncolumns++;
+	}
+
+	public void incrementNrows() {
+		nrows++;
+		yllCorner = yllCorner - cellSize;
+	}
+
+	public int getNcolumns() {
+		return ncolumns;
+	}
+
+	public int getNrows() {
+		return nrows;
+	}
+
+	public double getXllCorner() {
+		return xllCorner;
+	}
+
+	public double getYllCorner() {
+		return yllCorner;
+	}
+
+	public int getCellSize() {
+		return cellSize;
+	}
+
+	public String getNodataString() {
+		return nodataString;
+	}
+
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		result = prime * result + cellSize;
+		result = prime * result + ncolumns;
+		result = prime * result
+				+ ((nodataString == null) ? 0 : nodataString.hashCode());
+		result = prime * result + nrows;
+		long temp;
+		temp = Double.doubleToLongBits(xllCorner);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		temp = Double.doubleToLongBits(yllCorner);
+		result = prime * result + (int) (temp ^ (temp >>> 32));
+		return result;
+	}
+
+	@Override
+	public boolean equals(Object obj) {
+		if (this == obj)
+			return true;
+		if (obj == null)
+			return false;
+		if (!(obj instanceof RasterHeaderDetails))
+			return false;
+		RasterHeaderDetails other = (RasterHeaderDetails) obj;
+		if (cellSize != other.cellSize)
+			return false;
+		if (ncolumns != other.ncolumns)
+			return false;
+		if (nodataString == null) {
+			if (other.nodataString != null)
+				return false;
+		} else if (!nodataString.equals(other.nodataString))
+			return false;
+		if (nrows != other.nrows)
+			return false;
+		if (Double.doubleToLongBits(xllCorner) != Double
+				.doubleToLongBits(other.xllCorner))
+			return false;
+		if (Double.doubleToLongBits(yllCorner) != Double
+				.doubleToLongBits(other.yllCorner))
+			return false;
+		return true;
+	}
+	
+	@Override
+	public String toString() {
+		return "ncolumns=" + ncolumns + 
+				", nrows=" + nrows + 
+				", xllCorner=" + xllCorner + 
+				", yllCorner=" + yllCorner + 
+				", cellSize" + cellSize + 
+				", nodataString=" + nodataString;
+	}
+
+}
diff --git a/src/ac/sac/raster/RasterItem.java b/src/ac/sac/raster/RasterItem.java
new file mode 100755
index 00000000..42f30e02
--- /dev/null
+++ b/src/ac/sac/raster/RasterItem.java
@@ -0,0 +1,7 @@
+package ac.sac.raster;
+
+/*
+ * Marker interface for raster data.
+ * 
+ */
+public interface RasterItem {}
diff --git a/src/ac/sac/raster/RasterKey.java b/src/ac/sac/raster/RasterKey.java
new file mode 100755
index 00000000..2fc3a4e7
--- /dev/null
+++ b/src/ac/sac/raster/RasterKey.java
@@ -0,0 +1,69 @@
+package ac.sac.raster;
+
+import java.io.Serializable;
+
+public class RasterKey implements Serializable {
+
+	private static final long serialVersionUID = 7798904172536867819L;
+	private int col;
+	private int row;
+	
+	public RasterKey(RasterKey key) {
+		col = key.getCol();
+		row = key.getRow();
+	}
+	
+	public RasterKey(int col, int row) {
+		this.col = col;
+		this.row = row;
+	}
+	
+	public RasterKey createShifted(int colShift, int rowShift) {
+		if (colShift == 0 && rowShift == 0)
+			return this;
+		else
+			return new RasterKey(col + colShift, row + rowShift);
+	}
+	
+	public int getRow() {
+		return row;
+	}
+
+	public int getCol() {
+		return col;
+	}
+	
+	public double getDistanceTo(RasterKey x) {
+		double d = Math.sqrt(Math.pow(col - x.col,2) + Math.pow(row - x.row,2));
+		return d;
+	}
+	
+	@Override
+	public int hashCode() {
+		final int prime = 31;
+		int result = 1;
+		result = prime * result + col;
+		result = prime * result + row;
+		return result;
+	}
+	
+	@Override
+	public boolean equals(Object obj) {
+		if (this == obj)
+			return true;
+		if (obj == null)
+			return false;
+		if (!(obj instanceof RasterKey))
+			return false;
+		RasterKey other = (RasterKey) obj;
+		if (col != other.col)
+			return false;
+		if (row != other.row)
+			return false;
+		return true;
+	}	
+	
+	public String toString() {
+		return "row:" + row + ", col:" + col;
+	}
+}
diff --git a/src/ac/sac/raster/RasterSet.java b/src/ac/sac/raster/RasterSet.java
new file mode 100755
index 00000000..0114de81
--- /dev/null
+++ b/src/ac/sac/raster/RasterSet.java
@@ -0,0 +1,113 @@
+package ac.sac.raster;
+
+import java.util.ArrayList;
+import java.util.Collection;
+import java.util.HashMap;
+
+import ac.ed.lurg.utils.LogWriter;
+
+/* Class which holds raster data.  The generics used to defines the type of raster data held.
+ */
+public class RasterSet<D extends RasterItem> extends HashMap<RasterKey, D> {
+	
+	private static final long serialVersionUID = 4180188703734898215L;
+	private RasterHeaderDetails header;
+	
+	public RasterSet() {
+		header = new RasterHeaderDetails();
+	}
+	
+	public RasterSet(RasterHeaderDetails header) {
+		this.header = header;
+	}
+		
+	/** Returns the raster item for the row and column. If an source header is provided these are assumed to be that format system */
+	public Collection<D> getCollection(int col, int row, RasterHeaderDetails source) {
+		
+		Collection<D> dataPoints = new ArrayList<D>();
+		int c, r;
+		
+		if (source!=null && !source.equals(header)) {
+			int xOffset = (int) Math.round((source.getXllCorner() - header.getXllCorner())/header.getCellSize());
+			int yOffset = (int) Math.round((-source.getYllCorner() + header.getYllCorner())/header.getCellSize());
+
+			if (source.getCellSize() % header.getCellSize() != 0)
+				throw new RuntimeException("Cell sizes must be the same or souce an integer multiple of destination");
+			
+			int cellRatio = source.getCellSize() / header.getCellSize();
+			
+			for (int x = 0; x<cellRatio; x++) {
+				for (int y = 0; y<cellRatio; y++) {
+					c = (col)*cellRatio + x + xOffset;
+					r = (row)*cellRatio + y + yOffset;
+
+					D d = get(c, r);
+					if (d != null)
+						dataPoints.add(d);
+				}
+			}
+		}
+		else {
+			D d = get(col, row);
+			if (d != null)
+				dataPoints.add(d);
+		}
+		
+		return dataPoints;
+	}
+
+	/** Method to get raster item for coordinates */
+	public D getFromCoordinates(double source_x, double source_y) {
+		int col = (int) Math.round((source_x - header.getXllCorner())/header.getCellSize());
+		
+		// no idea why you need to turn this upside down by taking header.getNrows() and then subtracting, but it works
+		int row = header.getNrows() - (int) Math.round((source_y - header.getYllCorner())/header.getCellSize());
+		
+		if (row < 0 || col < 0) {
+			LogWriter.println("Got negative values");
+		}
+		return get(col, row);
+	}
+	
+	/** Method to really get the data, or create it. Assumes the col and row are in the internal header format */
+	public D get(int col, int row) {
+		if (header.getNcolumns() < col+1)
+			header.incrementNcolumns();
+		if (header.getNrows() < row+1) {
+			header.incrementNrows();
+		}
+		
+		RasterKey key = new RasterKey(col, row);
+		D data = get(key);
+		if (data == null) {
+			data = createRasterData();
+			if (data != null)
+				put(key, data);
+		}
+		return data;
+	}
+	
+	
+	/** Some classes do not create, those that do will need to override */
+	protected D createRasterData() {
+		return null;
+	}
+
+	public int getNcolumns() {
+		return header.getNcolumns();
+	}
+
+	public int getNrows() {
+		return header.getNrows();
+	}	
+	
+	public RasterHeaderDetails getHeaderDetails() {
+		return header;
+	}
+
+	/** Check passed in header is consistent with one from this set.
+	 *   Actually only checks grid size as can shift rasters*/
+	public boolean isConstistentWithHeader(RasterHeaderDetails h) {
+		return h.getCellSize() == header.getCellSize();
+	}
+}
-- 
GitLab