/////////////////////////////////////////////////////////////////////////////////////// /// \file demoinput.cpp /// \brief LPJ-GUESS input module for a toy data set (for demonstration purposes) /// /// /// \author Ben Smith /// $Date: 2013-11-22 10:56:18 +0100 (Fri, 22 Nov 2013) $ /// /////////////////////////////////////////////////////////////////////////////////////// #include "config.h" #include "demoinput.h" #include "driver.h" #include "outputchannel.h" #include REGISTER_INPUT_MODULE("demo", DemoInput) // Anonymous namespace for variables and functions with file scope namespace { // File names for temperature, precipitation, sunshine and soil code driver files xtring file_temp,file_prec,file_sun,file_soil; // LPJ soil code int soilcode; /// Interpolates monthly data to quasi-daily values. void interp_climate(double* mtemp, double* mprec, double* msun, double* mdtr, double* dtemp, double* dprec, double* dsun, double* ddtr) { interp_monthly_means_conserve(mtemp, dtemp, date.is_leap(date.year)); interp_monthly_totals_conserve(mprec, dprec, date.is_leap(date.year), 0); interp_monthly_means_conserve(msun, dsun, date.is_leap(date.year), 0, 100); interp_monthly_means_conserve(mdtr, ddtr, date.is_leap(date.year), 0); } } // namespace DemoInput::DemoInput() : nyear(1) { // Declare instruction file parameters declare_parameter("nyear", &nyear, 1, 10000, "Number of simulation years to run after spinup"); } bool DemoInput::read_from_file(Coord coord, xtring fname, const char* format, double monthly[12], bool soil /* = false */) { double dlon, dlat; int elev; FILE* in = fopen(fname, "r"); if (!in) { fail("readenv: could not open %s for input", (char*)fname); } bool foundgrid = false; while (!feof(in) && !foundgrid) { if (!soil) { readfor(in, format, &dlon, &dlat, &elev, monthly); } else { readfor(in, format, &dlon, &dlat, &soilcode); } foundgrid = equal(coord.lon, dlon) && equal(coord.lat, dlat); } fclose(in); if (!foundgrid) { dprintf("readenv: could not find record for (%g,%g) in %s", coord.lon, coord.lat, (char*)fname); } return foundgrid; } bool DemoInput::readenv(Coord coord, long& seed) { // Searches for environmental data in driver temperature, precipitation, // sunshine and soil code files for the grid cell whose coordinates are given by // 'coord'. Data are written to arrays mtemp, mprec, msun and the variable // soilcode, which are defined as global variables in this file // The temperature, precipitation and sunshine files (Cramer & Leemans, // unpublished) should be in ASCII text format and contain one-line records for the // climate (mean monthly temperature, mean monthly percentage sunshine or total // monthly precipitation) of a particular 0.5 x 0.5 degree grid cell. Elevation is // also included in each grid cell record. // The following sample record from the temperature file: // " -4400 8300 293-298-311-316-239-105 -7 26 -3 -91-184-239-277" // corresponds to the following data: // longitude 44 deg (-=W) // latitude 83 deg (+=N) // elevation 293 m // mean monthly temperatures (deg C) -29.8 (Jan), -31.1 (Feb), ..., -27.7 (Dec) // The following sample record from the precipitation file: // " 12750 -200 223 190 165 168 239 415 465 486 339 218 162 149 180" // corresponds to the following data: // longitude 127.5 deg (+=E) // latitude 20 deg (-=S) // elevation 223 m // monthly precipitation sum (mm) 190 (Jan), 165 (Feb), ..., 180 (Dec) // The following sample record from the sunshine file: // " 2600 7000 293 0 20 38 37 31 28 28 25 21 17 9 7" // corresponds to the following data: // longitude 26 deg (+=E) // latitude 70 deg (+=N) // elevation 293 m // monthly mean %age of full sunshine 0 (Jan), 20 (Feb), ..., 7 (Dec) // The LPJ soil code file is in ASCII text format and contains one-line records for // each grid cell in the form: // // where = longitude as a floating point number (-=W, +=E) // = latitude as a floating point number (-=S, +=N) // = integer in the range 0 (no soil) to 9 (see function // soilparameters in driver module) // The fields in each record are separated by spaces double mtemp[12]; // monthly mean temperature (deg C) double mprec[12]; // monthly precipitation sum (mm) double msun[12]; // monthly mean percentage sunshine values double mwet[12]={31,28,31,30,31,30,31,31,30,31,30,31}; // number of rain days per month double mdtr[12]; // monthly mean diurnal temperature range (oC) for(int m=0; m<12; m++) { mdtr[m] = 0.; if (ifbvoc) { dprintf("WARNING: No data available for dtr in sample data set!\nNo daytime temperature correction for BVOC calculations applied."); } } bool gridfound = read_from_file(coord, file_temp, "f6.2,f5.2,i4,12f4.1", mtemp); if(gridfound) gridfound = read_from_file(coord, file_prec, "f6.2,f5.2,i4,12f4", mprec); if(gridfound) gridfound = read_from_file(coord, file_sun, "f6.2,f5.2,i4,12f3", msun); if(gridfound) gridfound = read_from_file(coord, file_soil, "f,f,i", msun, true); // msun is not used here: just dummy if(gridfound) { // Interpolate monthly values for environmental drivers to daily values // (relevant daily values will be sent to the framework each simulation // day in function getclimate, below) interp_climate(mtemp, mprec, msun, mdtr, dtemp, dprec, dsun, ddtr); // Recalculate precipitation values using weather generator // (from Dieter Gerten 021121) prdaily(mprec, dprec, mwet, seed); } return gridfound; } void DemoInput::init() { // DESCRIPTION // Initialises input (e.g. opening files), and reads in the gridlist // // Reads list of grid cells and (optional) description text from grid list file // This file should consist of any number of one-line records in the format: // [] double dlon,dlat; bool eof=false; xtring descrip; // Read list of grid coordinates and store in global Coord object 'gridlist' // Retrieve name of grid list file as read from ins file xtring file_gridlist=param["file_gridlist"].str; FILE* in_grid=fopen(file_gridlist,"r"); if (!in_grid) fail("initio: could not open %s for input",(char*)file_gridlist); gridlist.killall(); first_call = true; while (!eof) { // Read next record in file eof=!readfor(in_grid,"f,f,a#",&dlon,&dlat,&descrip); if (!eof && !(dlon==0.0 && dlat==0.0)) { // ignore blank lines at end (if any) Coord& c=gridlist.createobj(); // add new coordinate to grid list c.lon=dlon; c.lat=dlat; c.descrip=descrip; } } fclose(in_grid); // Retrieve specified CO2 value as read from ins file co2=param["co2"].num; // Retrieve specified N value as read from ins file ndep=param["ndep"].num; // Open landcover files landcover_input.init(); // Open management files management_input.init(); // Retrieve input file names as read from ins file file_temp=param["file_temp"].str; file_prec=param["file_prec"].str; file_sun=param["file_sun"].str; file_soil=param["file_soil"].str; // Set timers tprogress.init(); tmute.init(); tprogress.settimer(); tmute.settimer(MUTESEC); } bool DemoInput::getgridcell(Gridcell& gridcell) { // See base class for documentation about this function's responsibilities // Select coordinates for next grid cell in linked list bool gridfound = false; bool LUerror = false; // Make sure we use the first gridcell in the first call to this function, // and then step through the gridlist in subsequent calls. if (first_call) { gridlist.firstobj(); // Note that first_call is static, so this assignment is remembered // across function calls. first_call = false; } else gridlist.nextobj(); if (gridlist.isobj) { while(!gridfound) { // Retrieve coordinate of next grid cell from linked list Coord& c = gridlist.getobj(); // Load environmental data for this grid cell from files if(run_landcover) { LUerror = landcover_input.loadlandcover(gridlist.getobj().lon, gridlist.getobj().lat); if(!LUerror) LUerror = management_input.loadmanagement(gridlist.getobj().lon, gridlist.getobj().lat); } if (!LUerror) { gridfound = readenv(c, gridcell.seed); } else { gridlist.nextobj(); if(!gridlist.isobj) return false; } } dprintf("\nCommencing simulation for stand at (%g,%g)",gridlist.getobj().lon, gridlist.getobj().lat); if (gridlist.getobj().descrip!="") dprintf(" (%s)\n\n", (char*)gridlist.getobj().descrip); else dprintf("\n\n"); // Tell framework the coordinates of this grid cell gridcell.set_coordinates(gridlist.getobj().lon, gridlist.getobj().lat); // The insolation data will be sent (in function getclimate, below) // as percentage sunshine gridcell.climate.instype=SUNSHINE; // Tell framework the soil type of this grid cell soilparameters(gridcell.soiltype,soilcode); // For Windows shell - clear graphical output // (ignored on other platforms) clear_all_graphs(); return true; // simulate this stand } return false; // no more stands } void DemoInput::getlandcover(Gridcell& gridcell) { landcover_input.getlandcover(gridcell); landcover_input.get_land_transitions(gridcell); } bool DemoInput::getclimate(Gridcell& gridcell) { // See base class for documentation about this function's responsibilities double progress; Climate& climate = gridcell.climate; // Send environmental values for today to framework climate.dndep = ndep / (365.0 * 10000.0); climate.co2 = co2; climate.temp = dtemp[date.day]; climate.prec = dprec[date.day]; climate.insol = dsun[date.day]; // bvoc climate.dtr=ddtr[date.day]; // First day of year only ... if (date.day == 0) { // Return false if last year was the last for the simulation if (date.year==nyear_spinup+nyear) return false; // Progress report to user and update timer if (tmute.getprogress()>=1.0) { progress=(double)(gridlist.getobj().id*(nyear_spinup+nyear) +date.year)/(double)(gridlist.nobj*(nyear_spinup+nyear)); tprogress.setprogress(progress); dprintf("%3d%% complete, %s elapsed, %s remaining\n",(int)(progress*100.0), tprogress.elapsed.str,tprogress.remaining.str); tmute.settimer(MUTESEC); } } return true; } DemoInput::~DemoInput() { // Performs memory deallocation, closing of files or other "cleanup" functions. }