Skip to content
Snippets Groups Projects
Commit 8456f44a authored by Bart Arendarczyk's avatar Bart Arendarczyk
Browse files

Binary file reader for carbon and wood yield data.

parent 363fc645
Branches
Tags
No related merge requests found
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
......@@ -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";
......
......
......@@ -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;
......
......
......@@ -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_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_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_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_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));
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);
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.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.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.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.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.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);
......
......
......@@ -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) {
......@@ -33,7 +39,7 @@ 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)) {
......
......
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);
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please to comment