diff --git a/config.properties b/config.properties new file mode 100644 index 0000000000000000000000000000000000000000..8ac5a792d734f459a81d0ff8d0f4455672110a2f --- /dev/null +++ b/config.properties @@ -0,0 +1,38 @@ +BASE_DIR=C:/Users/s1924442/git/plumv2 + +YIELD_DIR=C:/Users/s1924442/Documents/PhD/LURG/rcp60 +YIELD_DIR_TOP=rcp60 + +OUTPUT_DIR = C:/Users/s1924442/Documents/PhD/PLUM outputs/forest_calib + +WOOD_AND_CARBON_DIR = C:/Users/s1924442/Documents/PhD/Carbon and wood yields/forestry_dummy_4 + +BASE_YEAR=2010 +START_TIMESTEP=0 +TIMESTEP_SIZE=1 +END_TIMESTEP=10 + +IS_CALIBRATION_RUN=true +END_FIRST_STAGE_CALIBRATION=0 + +GENERATE_NEW_YIELD_CLUSTERS=false + +YIELD_FILENAME=yield.out + +DEBUG_LIMIT_COUNTRIES=true +DEBUG_COUNTRY_NAME=Central America +GAMS_COUNTRY_TO_SAVE=Central America + +INIT_WOOD_PRICE=0.4 +INIT_WOOD_STOCK=1000 +INIT_CARBON_PRICE=0.0 +INFRASTRUCTURE_EXPANSION_COST=0.0 +MANAGED_FOREST_INCREASE_COST=0.2 +MANAGED_FOREST_DECREASE_COST=0.2 +FOREST_MANAGEMENT_COST=0.25 +FOREST_MAX_CHANGE=0.02 +LAND_CHANGE_COST=0.5 +NATURAL_AREA_VALUE=0.0 +VEGETATION_CLEARING_COST=0.05 +IS_FORESTRY_ON=true +IS_CARBON_ON=true diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java index 17c3ee4d873f4f472bd8cdb5d28619c437a6527d..dcc2bd409a9b085a49b8656a26fa1d3163be6b55 100755 --- a/src/ac/ed/lurg/ModelConfig.java +++ b/src/ac/ed/lurg/ModelConfig.java @@ -250,8 +250,8 @@ public class ModelConfig { public static final String FOREST_DIR = SPATIAL_DATA_DIR + File.separator + "forestry"; public static final String WOOD_AND_CARBON_DIR = getProperty("WOOD_AND_CARBON_DIR"); //public static final String CARBON_FLUX_FILE = FOREST_DIR + File.separator + "carbon_flux_"; - public static final String WOOD_YIELD_AGRI_FILENAME = "wood_yield_to_agri.out"; - public static final String WOOD_YIELD_FRST_FILENAME = "wood_yield_to_frst.out"; + public static final String WOOD_YIELD_AGRI_FILENAME = "wood_yield_to_agri.dat"; + public static final String WOOD_YIELD_FRST_FILENAME = "wood_yield_to_frst.dat"; public static final String CARBON_LUC_FILENAME = "carbon_luc.out"; public static final String CARBON_NEE_FILENAME = "carbon_nee.out"; public static final String NATURAL_FOREST_GROWTH_FILENAME = "growth_natural.out"; diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java index 6eca89f0ae4f023c26e5f8d7594cee00b2c09e29..321964bdeb0d8420e6e034aaf7395ebdb75c683a 100644 --- a/src/ac/ed/lurg/ModelMain.java +++ b/src/ac/ed/lurg/ModelMain.java @@ -28,6 +28,7 @@ import ac.ed.lurg.demand.CalorieManager; import ac.ed.lurg.demand.DemandManagerFromFile; import ac.ed.lurg.demand.DemandManagerSSP; import ac.ed.lurg.demand.ElasticDemandManager; +import ac.ed.lurg.forestry.WoodYieldItem; import ac.ed.lurg.forestry.WoodYieldRasterSet; import ac.ed.lurg.landuse.ConversionCostReader; import ac.ed.lurg.landuse.CropUsageData; @@ -62,6 +63,8 @@ import ac.ed.lurg.utils.FileWriterHelper; import ac.ed.lurg.utils.LogWriter; import ac.ed.lurg.yield.LPJYieldResponseMapReader; import ac.ed.lurg.yield.YieldRaster; +import ac.sac.raster.AbstractBinaryRasterReader; +import ac.sac.raster.AbstractTabularRasterReader; import ac.sac.raster.IntegerRasterItem; import ac.sac.raster.IntegerRasterReader; import ac.sac.raster.InterpolatingRasterSet; diff --git a/src/ac/ed/lurg/carbon/CarbonFluxReader.java b/src/ac/ed/lurg/carbon/CarbonFluxReader.java index 3c4d8c7ea39585d93de2c3709ca5757149e3c6f6..8fc5350c8f6ade46d7f3ffc060f4487f5d1ef422 100644 --- a/src/ac/ed/lurg/carbon/CarbonFluxReader.java +++ b/src/ac/ed/lurg/carbon/CarbonFluxReader.java @@ -8,6 +8,7 @@ import ac.ed.lurg.Timestep; import ac.ed.lurg.landuse.LandUseItem; import ac.ed.lurg.landuse.LccKey; import ac.ed.lurg.types.LandCoverType; +import ac.sac.raster.AbstractBinaryRasterReader; import ac.sac.raster.AbstractTabularRasterReader; import ac.sac.raster.RasterHeaderDetails; import ac.sac.raster.RasterKey; @@ -18,45 +19,51 @@ public class CarbonFluxReader { private static final String DELIMITER = "[ |\t]+"; private static final double CONVERSION_FACTOR = 10.0; // convert kgC/m2 to tC/ha private RasterHeaderDetails rasterProj; + private String[] header = new String[MIN_COLS]; public CarbonFluxReader(RasterHeaderDetails rasterProj) { this.rasterProj = rasterProj; + header[0] = "Lon"; + header[1] = "Lat"; + for (int i = 0; i < MIN_COLS - 1; i++) { + header[i] = Integer.toString(i); + } } public CarbonFluxRasterSet getCarbonFluxes(RasterSet<LandUseItem> landUseRaster, Timestep timestep) { CarbonFluxRasterSet cFluxData = new CarbonFluxRasterSet(rasterProj); - getCarbonFluxesConversion("cflux_conv_ntrl_to_past.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.PASTURE)); - getCarbonFluxesConversion("cflux_conv_ntrl_to_forC.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.CARBON_FOREST)); - getCarbonFluxesConversion("cflux_conv_ntrl_to_forT.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.TIMBER_FOREST)); - getCarbonFluxesConversion("cflux_conv_ntrl_to_crop.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.CROPLAND)); + getCarbonFluxesConversion("cflux_conv_ntrl_to_past.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.PASTURE)); + getCarbonFluxesConversion("cflux_conv_ntrl_to_forC.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.CARBON_FOREST)); + getCarbonFluxesConversion("cflux_conv_ntrl_to_forT.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.TIMBER_FOREST)); + getCarbonFluxesConversion("cflux_conv_ntrl_to_crop.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.NATURAL, LandCoverType.CROPLAND)); - getCarbonFluxesConversion("cflux_conv_past_to_ntrl.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.NATURAL)); - getCarbonFluxesConversion("cflux_conv_past_to_forC.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.CARBON_FOREST)); - getCarbonFluxesConversion("cflux_conv_past_to_forT.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.TIMBER_FOREST)); - getCarbonFluxesConversion("cflux_conv_past_to_crop.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.CROPLAND)); + getCarbonFluxesConversion("cflux_conv_past_to_ntrl.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.NATURAL)); + getCarbonFluxesConversion("cflux_conv_past_to_forC.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.CARBON_FOREST)); + getCarbonFluxesConversion("cflux_conv_past_to_forT.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.TIMBER_FOREST)); + getCarbonFluxesConversion("cflux_conv_past_to_crop.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.PASTURE, LandCoverType.CROPLAND)); - getCarbonFluxesConversion("cflux_conv_forC_to_ntrl.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.NATURAL)); - getCarbonFluxesConversion("cflux_conv_forC_to_past.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.PASTURE)); - getCarbonFluxesConversion("cflux_conv_forC_to_forT.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.TIMBER_FOREST)); - getCarbonFluxesConversion("cflux_conv_forC_to_crop.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.CROPLAND)); + getCarbonFluxesConversion("cflux_conv_forC_to_ntrl.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.NATURAL)); + getCarbonFluxesConversion("cflux_conv_forC_to_past.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.PASTURE)); + getCarbonFluxesConversion("cflux_conv_forC_to_forT.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.TIMBER_FOREST)); + getCarbonFluxesConversion("cflux_conv_forC_to_crop.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CARBON_FOREST, LandCoverType.CROPLAND)); - getCarbonFluxesConversion("cflux_conv_forT_to_ntrl.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.NATURAL)); - getCarbonFluxesConversion("cflux_conv_forT_to_past.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.PASTURE)); - getCarbonFluxesConversion("cflux_conv_forT_to_forC.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CARBON_FOREST)); - getCarbonFluxesConversion("cflux_conv_forT_to_forT.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST)); - getCarbonFluxesConversion("cflux_conv_forT_to_crop.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CROPLAND)); + getCarbonFluxesConversion("cflux_conv_forT_to_ntrl.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.NATURAL)); + getCarbonFluxesConversion("cflux_conv_forT_to_past.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.PASTURE)); + getCarbonFluxesConversion("cflux_conv_forT_to_forC.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CARBON_FOREST)); + getCarbonFluxesConversion("cflux_conv_forT_to_forT.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.TIMBER_FOREST)); + getCarbonFluxesConversion("cflux_conv_forT_to_crop.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.TIMBER_FOREST, LandCoverType.CROPLAND)); - getCarbonFluxesConversion("cflux_conv_crop_to_ntrl.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.NATURAL)); - getCarbonFluxesConversion("cflux_conv_crop_to_past.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.PASTURE)); - getCarbonFluxesConversion("cflux_conv_crop_to_forC.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST)); - getCarbonFluxesConversion("cflux_conv_crop_to_forT.out", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST)); + getCarbonFluxesConversion("cflux_conv_crop_to_ntrl.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.NATURAL)); + getCarbonFluxesConversion("cflux_conv_crop_to_past.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.PASTURE)); + getCarbonFluxesConversion("cflux_conv_crop_to_forC.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST)); + getCarbonFluxesConversion("cflux_conv_crop_to_forT.dat", cFluxData, landUseRaster, timestep, new LccKey(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST)); - getCarbonFluxesEcosystem("cflux_sts_ntrl.out", cFluxData, landUseRaster, timestep, LandCoverType.NATURAL); - getCarbonFluxesEcosystem("cflux_sts_past.out", cFluxData, landUseRaster, timestep, LandCoverType.PASTURE); - getCarbonFluxesEcosystem("cflux_sts_forC.out", cFluxData, landUseRaster, timestep, LandCoverType.CARBON_FOREST); - getCarbonFluxesEcosystem("cflux_sts_forT.out", cFluxData, landUseRaster, timestep, LandCoverType.TIMBER_FOREST); - getCarbonFluxesEcosystem("cflux_sts_crop.out", cFluxData, landUseRaster, timestep, LandCoverType.CROPLAND); + getCarbonFluxesEcosystem("cflux_sts_ntrl.dat", cFluxData, landUseRaster, timestep, LandCoverType.NATURAL); + getCarbonFluxesEcosystem("cflux_sts_past.dat", cFluxData, landUseRaster, timestep, LandCoverType.PASTURE); + getCarbonFluxesEcosystem("cflux_sts_forC.dat", cFluxData, landUseRaster, timestep, LandCoverType.CARBON_FOREST); + getCarbonFluxesEcosystem("cflux_sts_forT.dat", cFluxData, landUseRaster, timestep, LandCoverType.TIMBER_FOREST); + getCarbonFluxesEcosystem("cflux_sts_crop.dat", cFluxData, landUseRaster, timestep, LandCoverType.CROPLAND); return cFluxData; @@ -64,7 +71,7 @@ public class CarbonFluxReader { public void getCarbonFluxesConversion(String filename, CarbonFluxRasterSet cFluxData, RasterSet<LandUseItem> landUseRaster, Timestep timestep, LccKey lccKey) { - AbstractTabularRasterReader<CarbonFluxItem> cFluxReader = new AbstractTabularRasterReader<CarbonFluxItem>(DELIMITER, MIN_COLS, cFluxData) { + AbstractBinaryRasterReader<CarbonFluxItem> cFluxReader = new AbstractBinaryRasterReader<CarbonFluxItem>(header, MIN_COLS - 2, cFluxData) { protected void setData(RasterKey key, CarbonFluxItem item, Map<String, Double> rowValues) { Double[] fluxes = getArrayFromRowValues(rowValues); @@ -77,7 +84,7 @@ public class CarbonFluxReader { public void getCarbonFluxesEcosystem(String filename, CarbonFluxRasterSet cFluxData, RasterSet<LandUseItem> landUseRaster, Timestep timestep, LandCoverType lcType) { - AbstractTabularRasterReader<CarbonFluxItem> cFluxReader = new AbstractTabularRasterReader<CarbonFluxItem>(DELIMITER, MIN_COLS, cFluxData) { + AbstractBinaryRasterReader<CarbonFluxItem> cFluxReader = new AbstractBinaryRasterReader<CarbonFluxItem>(header, MIN_COLS - 2, cFluxData) { protected void setData(RasterKey key, CarbonFluxItem item, Map<String, Double> rowValues) { Double[] fluxes = getArrayFromRowValues(rowValues); diff --git a/src/ac/ed/lurg/forestry/WoodYieldReader.java b/src/ac/ed/lurg/forestry/WoodYieldReader.java index 467b23f600af6f5d0cd2099bfa6fe512202dc91e..a4abac72e476b1a05d90022b69a6d5f731d37608 100644 --- a/src/ac/ed/lurg/forestry/WoodYieldReader.java +++ b/src/ac/ed/lurg/forestry/WoodYieldReader.java @@ -9,19 +9,25 @@ import ac.ed.lurg.Timestep; import ac.ed.lurg.landuse.LandUseItem; import ac.ed.lurg.landuse.LccKey; import ac.ed.lurg.types.LandCoverType; +import ac.sac.raster.AbstractBinaryRasterReader; import ac.sac.raster.AbstractTabularRasterReader; import ac.sac.raster.RasterHeaderDetails; import ac.sac.raster.RasterKey; import ac.sac.raster.RasterSet; public class WoodYieldReader { - private static final int MIN_COLS = 253; - private static final String DELIMITER = "[ |\t]+"; + private static final int MIN_COLS = 263; private static final double CONVERSION_FACTOR = 10.0; // convert kgC/m2 to tC/ha private RasterHeaderDetails rasterProj; + private String[] header = new String[MIN_COLS]; public WoodYieldReader(RasterHeaderDetails rasterProj) { this.rasterProj = rasterProj; + header[0] = "Lon"; + header[1] = "Lat"; + for (int i = 0; i < MIN_COLS - 1; i++) { + header[i] = Integer.toString(i); + } } public WoodYieldRasterSet getWoodYields(RasterSet<LandUseItem> landUseRaster, Timestep timestep, double woodPrice) { @@ -32,8 +38,8 @@ public class WoodYieldReader { } public void getWoodYieldsForest(WoodYieldRasterSet woodYieldData, RasterSet<LandUseItem> landUseRaster, Timestep timestep, double woodPrice) { - - AbstractTabularRasterReader<WoodYieldItem> yieldReader = new AbstractTabularRasterReader<WoodYieldItem>(DELIMITER, MIN_COLS, woodYieldData) { + + AbstractBinaryRasterReader<WoodYieldItem> yieldReader = new AbstractBinaryRasterReader<WoodYieldItem>(header, MIN_COLS - 2, woodYieldData) { protected void setData(RasterKey key, WoodYieldItem item, Map<String, Double> rowValues) { if (!landUseRaster.containsKey(key)) { @@ -61,7 +67,7 @@ public class WoodYieldReader { public void getWoodYieldsAgri(WoodYieldRasterSet woodYieldData, RasterSet<LandUseItem> landUseRaster, Timestep timestep, double woodPrice) { - AbstractTabularRasterReader<WoodYieldItem> yieldReader = new AbstractTabularRasterReader<WoodYieldItem>(DELIMITER, MIN_COLS, woodYieldData) { + AbstractBinaryRasterReader<WoodYieldItem> yieldReader = new AbstractBinaryRasterReader<WoodYieldItem>(header, MIN_COLS - 2, woodYieldData) { protected void setData(RasterKey key, WoodYieldItem item, Map<String, Double> rowValues) { if (!landUseRaster.containsKey(key)) { diff --git a/src/ac/sac/raster/AbstractBinaryRasterReader.java b/src/ac/sac/raster/AbstractBinaryRasterReader.java new file mode 100644 index 0000000000000000000000000000000000000000..46fbfb0eedc67efd580584dc6943e1b53598fd9d --- /dev/null +++ b/src/ac/sac/raster/AbstractBinaryRasterReader.java @@ -0,0 +1,108 @@ +package ac.sac.raster; + +import java.io.BufferedInputStream; +import java.io.File; +import java.io.FileInputStream; +import java.io.IOException; +import java.io.InputStream; +import java.nio.ByteBuffer; +import java.util.Arrays; +import java.util.HashMap; +import java.util.Map; + +import ac.ed.lurg.utils.LogWriter; + +public abstract class AbstractBinaryRasterReader<D extends RasterItem> { + + private static final int X_COL = 0; // longitude + private static final int Y_COL = 1; // latitude + private int NUM_DATA_COLS; // number of data columns + private String[] header; + protected RasterSet<D> dataset; + + public AbstractBinaryRasterReader(String[] header, int numDataCols, RasterSet<D> dataset) { + super(); + this.header = header; + NUM_DATA_COLS = numDataCols; + this.dataset = dataset; + } + + public RasterSet<D> getRasterDataFromFile(String filePath) { + long startTime = System.currentTimeMillis(); + + InputStream inputStream = null; + BufferedInputStream reader = null; + + try { + inputStream = new FileInputStream(filePath); + reader = new BufferedInputStream(inputStream); + + byte[] line = new byte[(NUM_DATA_COLS + 2) * Double.BYTES]; + + while (reader.read(line) > 0) { + Double[] dline = byteArrayToDouble(line); + + double x = dline[X_COL]; + double y = dline[Y_COL]; + + if (x > 180 || x < -180 || y > 90 || y < -90) { + int foo2 = 1; + } + + D item = dataset.getFromCoordinates(x, y); + RasterKey key = dataset.getKeyFromCoordinates(x, y); + + Map<String, Double> rowValues = new HashMap<String, Double>(NUM_DATA_COLS); + for (int i = 0; i < NUM_DATA_COLS; i++) { + rowValues.put(header[i], dline[i + 2]); + } + + setData(key, item, rowValues); + } + + } + catch (Exception e) { + LogWriter.printlnError("Problem reading data file " + filePath); + LogWriter.print(e); + throw new RuntimeException(e); + } + finally { + if (inputStream != null) { + try { + inputStream.close(); + } catch (IOException e) { + LogWriter.print(e); + } + } + if (reader != null) { + try { + reader.close(); + } catch (IOException e) { + LogWriter.print(e); + } + } + } + LogWriter.println("Reading " + filePath + ", took " + (System.currentTimeMillis() - startTime) + " ms"); + return dataset; + + } + + private Double[] byteArrayToDouble(byte[] bytes) { + int n = NUM_DATA_COLS + 2; + if (n * Double.BYTES != bytes.length) { + LogWriter.printlnError("AbstractBinaryRasterReader: byte array not of expected length"); + } + Double[] d = new Double[n]; + for (int i = 0; i < n; i++) { + ByteBuffer buff = ByteBuffer.allocate(Double.BYTES); + byte[] doubleBytes = Arrays.copyOfRange(bytes, i * Double.BYTES, i * Double.BYTES + Double.BYTES); + buff.put(doubleBytes); + buff.flip(); + d[i] = buff.getDouble(); + } + + return(d); + } + + abstract protected void setData(RasterKey key, D item, Map<String, Double> rowValues); +}