From cfed97e05998761729968ca8a91d137c4d7bae4f Mon Sep 17 00:00:00 2001 From: Ruth Rutter Date: Mon, 29 Apr 2024 14:34:37 -0700 Subject: [PATCH 1/3] Adds util function for last calendar year of input file --- include/TEMUtilityFunctions.h | 1 + src/TEMUtilityFunctions.cpp | 44 +++++++++++++++++++++++++++++++++++ 2 files changed, 45 insertions(+) diff --git a/include/TEMUtilityFunctions.h b/include/TEMUtilityFunctions.h index ea9073100..31aa45f84 100644 --- a/include/TEMUtilityFunctions.h +++ b/include/TEMUtilityFunctions.h @@ -145,6 +145,7 @@ namespace temutil { const int y, const int x); int get_timeseries_start_year(const std::string &filename); + int get_timeseries_end_year(const std::string &filename); // draft - reads all timesteps co2 data std::vector get_timeseries(const std::string &filename, diff --git a/src/TEMUtilityFunctions.cpp b/src/TEMUtilityFunctions.cpp index 3ba6d9b93..6bf24700f 100644 --- a/src/TEMUtilityFunctions.cpp +++ b/src/TEMUtilityFunctions.cpp @@ -659,6 +659,50 @@ namespace temutil { return start_year; } + /** Return the calendar end year of a timeseries netcdf file. + * Assumes that input files have complete years + * Assumes that time units are "days since..." + */ + int get_timeseries_end_year(const std::string& fname){ + + int start_year = get_timeseries_start_year(fname); + + int ncid; + temutil::nc( nc_open(fname.c_str(), NC_NOWRITE, &ncid) ); + + //Information about the time *dimension* + int timeD; + size_t timeD_len; + temutil::nc( nc_inq_dimid(ncid, "time", &timeD) ); + temutil::nc( nc_inq_dimlen(ncid, timeD, &timeD_len) ); + + //Information about the time *variable* + int timeV; + temutil::nc( nc_inq_varid(ncid, "time", &timeV) ); + + size_t start[3]; + start[0] = timeD_len-1; //Last entry only + start[1] = 0; + start[2] = 0; + + size_t count[3]; + count[0] = 1; + count[1] = 1; + count[2] = 1; + + int time_value; + temutil::nc( nc_get_vara_int(ncid, timeV, start, count, &time_value) ); + + //Round up, because the time value will be for the beginning + // of the last time step and so will not divide evenly by 365. + //Casting one of the values to a double to force use of the + // proper operator/ + int year_count = ceil(double(time_value) / DINY); + + int end_year = start_year + year_count; + return end_year; + } + /** rough draft - look up lon/lat in nc file from y,x coordinates. Assumes that the file has some coordinate dimensions... */ From 30fece937a583b8abba01f6bf687b8fcb6e2083a Mon Sep 17 00:00:00 2001 From: Ruth Rutter Date: Mon, 29 Apr 2024 15:25:16 -0700 Subject: [PATCH 2/3] Moves year range selection for EQ climate to config The calendar year values for historic climate averaging to produce EQ climate will now be specified in the configuration file. There is some basic checking to ensure that the values do not fall outside the historic input file range. --- config/config.js | 4 +++- include/Climate.h | 12 +++++++++-- include/ModelData.h | 4 ++++ src/Climate.cpp | 49 ++++++++++++++++++++++++++++++++++++--------- src/Cohort.cpp | 18 ++++++++++------- src/ModelData.cpp | 2 ++ 6 files changed, 69 insertions(+), 20 deletions(-) diff --git a/config/config.js b/config/config.js index 7fbd39726..151682ef9 100644 --- a/config/config.js +++ b/config/config.js @@ -121,7 +121,9 @@ }, "model_settings": { - "dynamic_lai": 1 // from model (1) or from input (0) + "dynamic_lai": 1, // from model (1) or from input (0) + "baseline_start": 1901, //start year for baseline EQ/SP climate + "baseline_end": 1931 //end year for baseline EQ/SP climate // //"dynamic_climate": 0, // //"varied_co2": 0, // //"fire_severity_as_input": 0, // fire sev. as input or ?? diff --git a/include/Climate.h b/include/Climate.h index 7a5a81866..279f41634 100644 --- a/include/Climate.h +++ b/include/Climate.h @@ -10,10 +10,16 @@ class Climate { Climate(); Climate(const std::string& fname, const std::string& co2fname, int y, int x); - // misc climate variables - // This value will change during a run when switching from + // Misc. climate variables + // The following two values determine the span of historic climate + // years used to produce a baseline average climate for EQ + // TODO: SP repeated climate is not currently affected by these. + int baseline_start; + int baseline_end; + // These values will change during a run when switching from // historic to projected climate data int tseries_start_year; + int tseries_end_year; // driving variables std::vector co2; @@ -75,6 +81,8 @@ class Climate { void load_proj_climate(const std::string&, int, int); void load_proj_co2(const std::string&); + void prep_avg_climate(); + private: void load_from_file(const std::string& fname, int y, int x); diff --git a/include/ModelData.h b/include/ModelData.h index 95e00898f..43e3a5319 100644 --- a/include/ModelData.h +++ b/include/ModelData.h @@ -86,6 +86,10 @@ class ModelData { bool nc_tr; bool nc_sc; int output_interval; //How many years to store for output + //The following two config values are temporarily stored in + // ModelData, to be transferred to Climate. + int baseline_start;//Start year for baseline EQ climate + int baseline_end;//End year for baseline EQ climate // Maps holding data about variables to be output at specific timesteps // C++11 would allow the use of unordered_maps, which have a faster diff --git a/src/Climate.cpp b/src/Climate.cpp index 6b240c367..a8164d531 100644 --- a/src/Climate.cpp +++ b/src/Climate.cpp @@ -405,6 +405,7 @@ void Climate::load_from_file(const std::string& fname, int y, int x) { nirr = temutil::get_timeseries(fname, "nirr", y, x); tseries_start_year = temutil::get_timeseries_start_year(fname); + tseries_end_year = temutil::get_timeseries_end_year(fname); }//End critical(load_climate) @@ -454,14 +455,15 @@ void Climate::load_from_file(const std::string& fname, int y, int x) { BOOST_LOG_SEV(glg, debug) << "par = [" << temutil::vec2csv(par) << "]"; // Create a simplified historic climate for EQ by averaging input data - // over a year range specified here. The year choices should be exposed - // in the config file eventually. TODO - if(fname.find("historic") != std::string::npos){ - avgX_tair = avg_over(tair, 1901, 1931); - avgX_prec = avg_over(prec, 1901, 1931); - avgX_nirr = avg_over(nirr, 1901, 1931); - avgX_vapo = avg_over(vapo, 1901, 1931); - } + // over a year range specified here. +//Commenting out but leaving temporarily for easier comparison with +// the currently in-progress daily data ingestion branch. 20240429 +// if(fname.find("historic") != std::string::npos){ +// avgX_tair = avg_over(tair, baseline_start, baseline_end); +// avgX_prec = avg_over(prec, 1901, 1931); +// avgX_nirr = avg_over(nirr, 1901, 1931); +// avgX_vapo = avg_over(vapo, 1901, 1931); +// } // Do we need simplified 'avgX_' values for par, and cld?? // ===> YES: the derived variables should probably be based off the avgX @@ -474,6 +476,35 @@ void Climate::load_from_file(const std::string& fname, int y, int x) { } +/** Prepares simplified average climate, intended for EQ stage +* +* Assumes that averaged climate is produced shortly after +* the tseries values have been loaded from the intended +* input file. +*/ +void Climate::prep_avg_climate(){ + BOOST_LOG_SEV(glg, debug) << "Climate baseline from config: " + << this->climate.baseline_start << ":" + << this->climate.baseline_end; + + if( baseline_start > tseries_end_year + || baseline_start < tseries_start_year + || baseline_end > tseries_end_year + || baseline_end < tseries_start_year){ + + BOOST_LOG_SEV(glg, fatal) << "Baseline year value exceeds range" + << " of input file.\n" << "Acceptable range: " + << tseries_start_year << ":" << tseries_end_year << std::endl; + exit(EXIT_FAILURE); + } + + avgX_tair = avg_over(tair, baseline_start, baseline_end); + avgX_prec = avg_over(prec, baseline_start, baseline_end); + avgX_nirr = avg_over(nirr, baseline_start, baseline_end); + avgX_vapo = avg_over(vapo, baseline_start, baseline_end); +} + + /** This loads data from a projected climate data file, overwriting any old climate data*/ void Climate::load_proj_climate(const std::string& fname, int y, int x){ BOOST_LOG_SEV(glg, note) << "Climate, loading projected data"; @@ -518,9 +549,7 @@ std::vector Climate::avg_over(const std::vector & var, const int s } BOOST_LOG_SEV(glg, debug) << "result = [" << temutil::vec2csv(result) << "]"; - return result; - } diff --git a/src/Cohort.cpp b/src/Cohort.cpp index ec8e4b0b9..48ab27266 100644 --- a/src/Cohort.cpp +++ b/src/Cohort.cpp @@ -63,13 +63,17 @@ Cohort::Cohort(int y, int x, ModelData* modeldatapointer): // might need to set the cd* and the ed* ??? - BOOST_LOG_SEV(glg, debug) << "Setup the NEW STYLE CLIMATE OBJECT ..."; - // FIX: Historic? Projected?? how to handle both?? - // Maybe: - //this->hist_climate = Climate(modeldatapointer->hist_climate, y, x); - //this->proj_climate = Climate(modeldatapointer->proj_climate, y, x); + BOOST_LOG_SEV(glg, debug) << "Construct Climate object ..."; + + //On construction we assume historic climate, which will be + // overwritten by projected climate later when necessary. this->climate = Climate(modeldatapointer->hist_climate_file, modeldatapointer->co2_file, y, x); - + + this->climate.baseline_start = modeldatapointer->baseline_start; + this->climate.baseline_end = modeldatapointer->baseline_end; + //Prepare averaged input set for EQ stage + this->climate.prep_avg_climate(); + // Build a mineral info object MineralInfo mineral_info = MineralInfo(modeldatapointer->soil_texture_file, y, x); @@ -853,6 +857,7 @@ void Cohort::updateMonthly_DIMveg(const int & currmind, const bool & dynamic_lai for (int ip=0; ip0.) { cd.m_veg.vegage[ip] = cd.yrsdist; + if (cd.m_veg.vegage[ip]<=0) { cd.m_vegd.foliagemx[ip] = 0.; } @@ -1518,4 +1523,3 @@ void Cohort::set_restartdata_from_state() { } } } - diff --git a/src/ModelData.cpp b/src/ModelData.cpp index 182a4b27b..24d2b6dbf 100644 --- a/src/ModelData.cpp +++ b/src/ModelData.cpp @@ -112,6 +112,8 @@ ModelData::ModelData(Json::Value controldata):force_cmt(-1) { caldata_tree_loc = controldata["calibration-IO"]["caldata_tree_loc"].asString(); dynamic_LAI = controldata["model_settings"]["dynamic_lai"].asInt(); // checked in Cohort::updateMonthly_DIMVeg + baseline_start = controldata["model_settings"]["baseline_start"].asInt(); + baseline_end = controldata["model_settings"]["baseline_end"].asInt(); // Unused (11/23/2015) //changeclimate = controldata["model_settings"]["dynamic_climate"].asInt(); From a74758a1b0b551a944d241bdc7e2a4797479f6d3 Mon Sep 17 00:00:00 2001 From: Ruth Rutter Date: Mon, 29 Apr 2024 15:54:07 -0700 Subject: [PATCH 3/3] Adds more input checking, fixes some details --- config/config.js | 4 ++-- src/Climate.cpp | 10 ++++++++-- 2 files changed, 10 insertions(+), 4 deletions(-) diff --git a/config/config.js b/config/config.js index 151682ef9..49fd6b912 100644 --- a/config/config.js +++ b/config/config.js @@ -122,8 +122,8 @@ "model_settings": { "dynamic_lai": 1, // from model (1) or from input (0) - "baseline_start": 1901, //start year for baseline EQ/SP climate - "baseline_end": 1931 //end year for baseline EQ/SP climate + "baseline_start": 1901, //start year for baseline EQ climate + "baseline_end": 1931 //end year for baseline EQ climate // //"dynamic_climate": 0, // //"varied_co2": 0, // //"fire_severity_as_input": 0, // fire sev. as input or ?? diff --git a/src/Climate.cpp b/src/Climate.cpp index a8164d531..14ca54bfd 100644 --- a/src/Climate.cpp +++ b/src/Climate.cpp @@ -484,8 +484,7 @@ void Climate::load_from_file(const std::string& fname, int y, int x) { */ void Climate::prep_avg_climate(){ BOOST_LOG_SEV(glg, debug) << "Climate baseline from config: " - << this->climate.baseline_start << ":" - << this->climate.baseline_end; + << baseline_start << ":" << baseline_end; if( baseline_start > tseries_end_year || baseline_start < tseries_start_year @@ -498,6 +497,13 @@ void Climate::prep_avg_climate(){ exit(EXIT_FAILURE); } + if(baseline_end < baseline_start){ + BOOST_LOG_SEV(glg, fatal) << "Baseline end year " << baseline_end + << " is less than baseline start year " + << baseline_start; + exit(EXIT_FAILURE); + } + avgX_tair = avg_over(tair, baseline_start, baseline_end); avgX_prec = avg_over(prec, baseline_start, baseline_end); avgX_nirr = avg_over(nirr, baseline_start, baseline_end);