From 538698379f3ec7d6474f197c7bebc2f1dacabdc0 Mon Sep 17 00:00:00 2001
From: Peter Alexander <peter@blackhillock.co.uk>
Date: Tue, 4 Aug 2015 12:09:08 +0100
Subject: [PATCH] Produce area output files for LPJ

---
 src/ac/ed/lurg/ModelMain.java                 |  33 +++---
 src/ac/ed/lurg/landuse/AreasItem.java         |   9 +-
 src/ac/ed/lurg/output/LpjgOutputer.java       | 103 ++++++++++++------
 .../lurg/yield/LPJYieldResponseMapReader.java |  15 ++-
 src/ac/ed/lurg/yield/YieldResponsesItem.java  |   9 ++
 5 files changed, 114 insertions(+), 55 deletions(-)

diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 1e54e573..cbd31705 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -24,6 +24,7 @@ import ac.ed.lurg.landuse.IrrigationCostItem;
 import ac.ed.lurg.landuse.IrrigiationCostReader;
 import ac.ed.lurg.landuse.LandCoverItem;
 import ac.ed.lurg.landuse.LandCoverReader;
+import ac.ed.lurg.output.LpjgOutputer;
 import ac.ed.lurg.types.CropToDoubleMap;
 import ac.ed.lurg.types.CropType;
 import ac.ed.lurg.types.FertiliserRate;
@@ -77,7 +78,8 @@ public class ModelMain {
 		//		YieldResponsesItem yresp = yieldSurfaces.getFromCoordinates(-50.0, -4.0);
 		//		LogWriter.printlnError("Test key: " + yresp.getYieldMax(CropType.CEREALS)  + ", " + yresp.getYieldFertOnly(CropType.CEREALS) + ", " + yresp.getYieldIrrigOnly(CropType.CEREALS));
 
-		LogWriter.println("Timestep " + timestep);
+		int year = ModelConfig.BASE_YEAR + timestep;
+		LogWriter.println("Timestep: " + timestep + ", year:" + year);
 
 		CropToDoubleMap totalQuantity = new CropToDoubleMap();
 		CropToDoubleMap totalWorldInputCost = new CropToDoubleMap();
@@ -133,12 +135,15 @@ public class ModelMain {
 		prevWorldInputCost = totalWorldInputCost.divideBy(totalQuantity);
 		
 		// output results
-		outputTimestepResults(timestep, globalIntensityRaster, globalCropAreaRaster, globalLocationIdRaster);
+		outputTimestepResults(year, globalIntensityRaster, globalCropAreaRaster, globalLocationIdRaster, yieldSurfaces);
 	}
 
-	private void outputTimestepResults(int timestep, RasterSet<IntensitiesItem> intensityRaster, RasterSet<AreasItem> cropAreaRaster, RasterSet<IntegerRasterItem> locationIdRaster) {
+	private void outputTimestepResults(int year, RasterSet<IntensitiesItem> intensityRaster, RasterSet<AreasItem> cropAreaRaster, RasterSet<IntegerRasterItem> locationIdRaster, YieldRaster yieldSurfaces) {
 		
-		new RasterOutputer<IntegerRasterItem>(locationIdRaster, "locId" + timestep) {
+		LpjgOutputer lpjOutputer = new LpjgOutputer(year, intensityRaster, cropAreaRaster, yieldSurfaces);
+		lpjOutputer.writeOutput();
+		
+		/*	new RasterOutputer<IntegerRasterItem>(locationIdRaster, "locId" + year) {
 			@Override
 			public Double getValue(RasterKey location) {
 				IntegerRasterItem locItem = results.get(location);
@@ -155,7 +160,7 @@ public class ModelMain {
 		}.writeOutput(ModelConfig.WRITE_JPEG_IMAGES);
 
 
-		new RasterOutputer<IntensitiesItem>(intensityRaster, "wheatIntensity" + timestep) {
+		new RasterOutputer<IntensitiesItem>(intensityRaster, "wheatIntensity" + year) {
 			@Override
 			public Double getValue(RasterKey location) {
 				IntensitiesItem intensityItem = results.get(location);
@@ -173,13 +178,13 @@ public class ModelMain {
 			}	
 		}.writeOutput(ModelConfig.WRITE_JPEG_IMAGES);
 
-		outputAreas(timestep, cropAreaRaster, CropType.WHEAT);
-		outputAreas(timestep, cropAreaRaster, CropType.MAIZE);
-		outputAreas(timestep, cropAreaRaster, CropType.PASTURE);
+		outputAreas(year, cropAreaRaster, CropType.MAIZE); */
+		outputAreas(year, cropAreaRaster, CropType.WHEAT);
+		outputAreas(year, cropAreaRaster, CropType.PASTURE); 
 	}
 
-	private void outputAreas(int timestep, RasterSet<AreasItem> cropAreaRaster, final CropType crop) {
-		new RasterOutputer<AreasItem>(cropAreaRaster, crop.getGamsName() + "Area" + timestep) {
+	private void outputAreas(int year, RasterSet<AreasItem> cropAreaRaster, final CropType crop) {
+		new RasterOutputer<AreasItem>(cropAreaRaster, crop.getGamsName() + "Area" + year) {
 			@Override
 			public Double getValue(RasterKey location) {
 				AreasItem area = results.get(location);
@@ -241,11 +246,11 @@ public class ModelMain {
 
 
 			// DEBUG code
-			if (!(country.getCountryName().equals("United States of Americaxx") || country.getCountryName().equals("Russian Federationxx") || country.getCountryName().equals("China")) ) { //|| country.getCountryName().equals("China")
-				continue;
-			}
+		//	if (!(country.getCountryName().equals("United States of Americaxx") || country.getCountryName().equals("Russian Federationxx") || country.getCountryName().equals("China")) ) { //|| country.getCountryName().equals("China")
+		//		continue;
+		//	}
 
-			if (demandManager.getPopulation(country, 2010) < 10 || countryExclusionList.contains(country.getCountryName())) {
+			if (demandManager.getPopulation(country, 2010) < 50 || countryExclusionList.contains(country.getCountryName())) {
 				LogWriter.printlnError("Skipping " + country);
 				continue;
 			}
diff --git a/src/ac/ed/lurg/landuse/AreasItem.java b/src/ac/ed/lurg/landuse/AreasItem.java
index e5b6f4d7..6a24ed18 100644
--- a/src/ac/ed/lurg/landuse/AreasItem.java
+++ b/src/ac/ed/lurg/landuse/AreasItem.java
@@ -54,6 +54,14 @@ public class AreasItem implements RasterItem {
 		return d == null ? 0.0 : d;
 	}
 
+	public double getCropFraction(CropType... cropsToFind) {
+		double totalFract = 0;
+		for (CropType crop : cropsToFind) {
+			totalFract += getCropFraction(crop);
+		}
+		return totalFract;
+	}
+
 	
 	public void setCropFraction(CropType c, double area) {
 		cropFractions.put(c, area);
@@ -93,5 +101,4 @@ public class AreasItem implements RasterItem {
 		
 		return changes;
 	}
-
 }
diff --git a/src/ac/ed/lurg/output/LpjgOutputer.java b/src/ac/ed/lurg/output/LpjgOutputer.java
index 5a3c8908..718ca692 100644
--- a/src/ac/ed/lurg/output/LpjgOutputer.java
+++ b/src/ac/ed/lurg/output/LpjgOutputer.java
@@ -9,6 +9,11 @@ import java.util.Map.Entry;
 import ac.ed.lurg.ModelConfig;
 import ac.ed.lurg.landuse.AreasItem;
 import ac.ed.lurg.landuse.IntensitiesItem;
+import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.utils.LogWriter;
+import ac.ed.lurg.yield.YieldRaster;
+import ac.ed.lurg.yield.YieldResponsesItem;
 import ac.sac.raster.RasterKey;
 import ac.sac.raster.RasterSet;
 
@@ -16,54 +21,84 @@ public class LpjgOutputer {
 
 	private RasterSet<IntensitiesItem> intensityRaster;
 	private RasterSet<AreasItem> cropAreaRaster;
-	private int timestep;
-	
-	public LpjgOutputer(int timestep, RasterSet<IntensitiesItem> intensityRaster, RasterSet<AreasItem> cropAreaRaster) {
+	private YieldRaster yieldSurfaces;
+	private int year;
+
+	public LpjgOutputer(int year, RasterSet<IntensitiesItem> intensityRaster, RasterSet<AreasItem> cropAreaRaster, YieldRaster yieldSurfaces) {
+		this.year = year;
 		this.intensityRaster = intensityRaster;
 		this.cropAreaRaster = cropAreaRaster;
+		this.yieldSurfaces = yieldSurfaces;
+	}
+
+	public void writeOutput() {
+		writeLandCoverAndCrop();
 	}
 	
-	
-	public void writeOutput(boolean writeImage) {
-	    BufferedWriter fileWriter = null;
+	public void writeLandCoverAndCrop() {
+		BufferedWriter landCoverWriter = null;
+		BufferedWriter cropFractWriter = null;
 
-	    try {
-    		String areaFileName = ModelConfig.OUTPUT_DIR + File.separator + "Areas" + timestep + ".txt";
-    		fileWriter = new BufferedWriter(new FileWriter(areaFileName,false));
-    
-	    	writeHeader(fileWriter);
+		try {
+			String landCoverFileName = ModelConfig.OUTPUT_DIR + File.separator + "LandCoverFract" + year + ".txt";
+			landCoverWriter = new BufferedWriter(new FileWriter(landCoverFileName, false));
+			landCoverWriter.write("Lon Lat CROPLAND PASTURE NATURAL");
+			landCoverWriter.newLine();
 
+			String cropFractionFileName = ModelConfig.OUTPUT_DIR + File.separator + "CropFract" + year + ".txt";
+			cropFractWriter = new BufferedWriter(new FileWriter(cropFractionFileName, false));
+			cropFractWriter.write("Lon Lat TeWWirr TeSWirr TeCoirr TrRiirr");
+			cropFractWriter.newLine();
+			
 			for (Entry<RasterKey, AreasItem> entry : cropAreaRaster.entrySet()) {
-					RasterKey location = entry.getKey();
-					AreasItem item = entry.getValue();
+				RasterKey key = entry.getKey();
+				AreasItem item = entry.getValue();
+				double lat = cropAreaRaster.getXCoordin(key);
+				double lon = cropAreaRaster.getYCoordin(key);
+
+				double area = cropAreaRaster.getAreaMha(key);
+				double crop = item.getLandCoverArea(LandCoverType.CROPLAND);
+				double pasture = item.getLandCoverArea(LandCoverType.PASTURE);
+				double forest = item.getLandCoverArea(LandCoverType.FOREST);
+				double otherNatural = item.getLandCoverArea(LandCoverType.OTHER_NATURAL);
+				landCoverWriter.write(String.format("%.2f %.2f %.4f %.4f %.4f", lat, lon, crop/area, pasture/area, (forest+otherNatural)/area));
+				landCoverWriter.newLine();
 
-		/*			Double d = getValue(location);
-					
-					if (d == null)
-						fileWriter.write(nullDataString + " ");
-					else
-						fileWriter.write(d + " ");*/
+				boolean isSpringWheat = false;
+				YieldResponsesItem yieldItem = yieldSurfaces.get(key);
+				if (yieldItem != null)
+					isSpringWheat = yieldSurfaces.get(key).isSpringWheat();
+				else 
+					LogWriter.printlnError("LpjgOutputer: Can't find YieldResponsesItem for " + key);
 				
-				fileWriter.newLine();
-	    	}
-	    }
-	    catch (IOException e) {
-	    	e.printStackTrace();
-	    }
-	    finally {
-			if (fileWriter != null) {
+				double lpjWinterWheatFrac = item.getCropFraction(CropType.OILCROPS, isSpringWheat ? null : CropType.WHEAT);
+				double lpjSpringWheatFrac = item.getCropFraction(CropType.SOYBEAN, CropType.PULSES, CropType.STARCHY_ROOTS, isSpringWheat ? CropType.WHEAT : null);
+				double lpjCornFrac = item.getCropFraction(CropType.MAIZE, CropType.SOYBEAN);
+				double lpjRiceFrac = item.getCropFraction(CropType.RICE);
+				cropFractWriter.write(String.format("%.2f %.2f %.4f %.4f %.4f %.4f", lat, lon, lpjWinterWheatFrac, lpjSpringWheatFrac, lpjCornFrac, lpjRiceFrac));
+				cropFractWriter.newLine();
+			}
+		}
+		catch (IOException e) {
+			e.printStackTrace();
+		}
+		finally {
+			if (landCoverWriter != null) {
 				try {
-					fileWriter.close();
+					landCoverWriter.close();
 				} 
 				catch (IOException e) {
 					e.printStackTrace();
 				}
 			}
-	    }
-	}
-
-	private void writeHeader(BufferedWriter outputFile) throws IOException {
-		outputFile.write("h1,h2,... ,hn");
-    	outputFile.newLine();
+			if (cropFractWriter != null) {
+				try {
+					cropFractWriter.close();
+				} 
+				catch (IOException e) {
+					e.printStackTrace();
+				}
+			}
+		}
 	}
 }
diff --git a/src/ac/ed/lurg/yield/LPJYieldResponseMapReader.java b/src/ac/ed/lurg/yield/LPJYieldResponseMapReader.java
index bcdea941..5a7e1573 100644
--- a/src/ac/ed/lurg/yield/LPJYieldResponseMapReader.java
+++ b/src/ac/ed/lurg/yield/LPJYieldResponseMapReader.java
@@ -112,8 +112,11 @@ public class LPJYieldResponseMapReader {
 	
 	public void setData(YieldType noIrrigYieldType, YieldType maxIrrigYieldType, YieldResponsesItem item, YieldRecord record) {
 		
-		item.setYield(noIrrigYieldType, CropType.WHEAT, Math.max(record.teSW, record.teWW));
-		item.setYield(maxIrrigYieldType, CropType.WHEAT, Math.max(record.teSWirr, record.teWWirr));
+		boolean isSpringWheat = (record.teSW > record.teWW);
+		item.setIsSpringWheat(isSpringWheat);
+		
+		item.setYield(noIrrigYieldType, CropType.WHEAT, isSpringWheat ? record.teSW : record.teWW);
+		item.setYield(maxIrrigYieldType, CropType.WHEAT, isSpringWheat ? record.teSWirr : record.teWWirr);
 		item.setYield(noIrrigYieldType, CropType.MAIZE, record.teCo);
 		item.setYield(maxIrrigYieldType, CropType.MAIZE, record.teCoirr);
 		item.setYield(noIrrigYieldType, CropType.RICE, record.trRi);
@@ -127,11 +130,11 @@ public class LPJYieldResponseMapReader {
 		item.setYield(maxIrrigYieldType, CropType.SOYBEAN, record.teSWirr);
 		item.setYield(noIrrigYieldType, CropType.PULSES, record.teSW);
 		item.setYield(maxIrrigYieldType, CropType.PULSES, record.teSWirr);
-		item.setYield(noIrrigYieldType, CropType.STARCHY_ROOTS, record.teSW * 1.1);
-		item.setYield(maxIrrigYieldType, CropType.STARCHY_ROOTS, record.teSWirr * 1.1);
+		item.setYield(noIrrigYieldType, CropType.STARCHY_ROOTS, record.teSW);
+		item.setYield(maxIrrigYieldType, CropType.STARCHY_ROOTS, record.teSWirr);
 		
-		item.setYield(noIrrigYieldType, CropType.PASTURE, Math.max(record.c3Pasture, record.c4Pasture));
-		item.setYield(maxIrrigYieldType, CropType.PASTURE, Math.max(record.c3Pasture, record.c4Pasture)*1.2);
+		item.setYield(noIrrigYieldType, CropType.PASTURE, Math.max(record.c3Pasture + record.c4Pasture, 0));
+		item.setYield(maxIrrigYieldType, CropType.PASTURE,  Math.max(record.c3Pasture + record.c4Pasture, 0)*1.2);
 	}
 	
 	class YieldRecord {
diff --git a/src/ac/ed/lurg/yield/YieldResponsesItem.java b/src/ac/ed/lurg/yield/YieldResponsesItem.java
index cc8f6195..723a12cd 100644
--- a/src/ac/ed/lurg/yield/YieldResponsesItem.java
+++ b/src/ac/ed/lurg/yield/YieldResponsesItem.java
@@ -10,6 +10,7 @@ import ac.sac.raster.RasterItem;
 
 public class YieldResponsesItem implements RasterItem {
 	Map<CropType, YieldResponse> yieldResponses = new HashMap<CropType, YieldResponse>();
+	private boolean isSpringWheat;  // used to indicate is wheat is spring or winter wheat
 	
 	public void setYield(YieldType type, CropType crop, double yield) {
 		getYieldResponseForCrop(crop).setYield(type, yield);
@@ -56,4 +57,12 @@ public class YieldResponsesItem implements RasterItem {
 	public double getIrrigParam(CropType crop) {
 		return getYieldResponseForCrop(crop).getIrrigParam();
 	}
+
+	public void setIsSpringWheat(boolean isSpringWheat) {
+		this.isSpringWheat = isSpringWheat;
+	}
+	
+	public boolean isSpringWheat() {
+		return isSpringWheat;
+	}
 }
\ No newline at end of file
-- 
GitLab