FF16 Rainfall API

R/RcppExports.R

@@ -3293,6 +3293,38 @@ FF16_Environment__ctor <- function(canopy_rescale_usually, soil_number_of_depths
     .Call('_plant_FF16_Environment__ctor', PACKAGE = 'plant', canopy_rescale_usually, soil_number_of_depths)
 }
 
+FF16_Environment__rainfall_init <- function(obj_, x, y) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_init', PACKAGE = 'plant', obj_, x, y))
+}
+
+FF16_Environment__rainfall_recompute <- function(obj_) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_recompute', PACKAGE = 'plant', obj_))
+}
+
+FF16_Environment__rainfall_eval <- function(obj_, u) {
+    .Call('_plant_FF16_Environment__rainfall_eval', PACKAGE = 'plant', obj_, u)
+}
+
+FF16_Environment__rainfall_eval_range <- function(obj_, u) {
+    .Call('_plant_FF16_Environment__rainfall_eval_range', PACKAGE = 'plant', obj_, u)
+}
+
+FF16_Environment__rainfall_add_point <- function(obj_, x, y) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_add_point', PACKAGE = 'plant', obj_, x, y))
+}
+
+FF16_Environment__rainfall_add_point_sorted <- function(obj_, x, y) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_add_point_sorted', PACKAGE = 'plant', obj_, x, y))
+}
+
+FF16_Environment__rainfall_add_points <- function(obj_, x, y) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_add_points', PACKAGE = 'plant', obj_, x, y))
+}
+
+FF16_Environment__rainfall_clear_points <- function(obj_) {
+    invisible(.Call('_plant_FF16_Environment__rainfall_clear_points', PACKAGE = 'plant', obj_))
+}
+
 FF16_Environment__canopy_openness <- function(obj_, height) {
     .Call('_plant_FF16_Environment__canopy_openness', PACKAGE = 'plant', obj_, height)
 }
@@ -3305,10 +3337,6 @@ FF16_Environment__set_fixed_environment <- function(obj_, value, height_max) {
     invisible(.Call('_plant_FF16_Environment__set_fixed_environment', PACKAGE = 'plant', obj_, value, height_max))
 }
 
-FF16_Environment__set_soil_infiltration_rate <- function(obj_, rate) {
-    invisible(.Call('_plant_FF16_Environment__set_soil_infiltration_rate', PACKAGE = 'plant', obj_, rate))
-}
-
 FF16_Environment__set_soil_water_state <- function(obj_, state) {
     invisible(.Call('_plant_FF16_Environment__set_soil_water_state', PACKAGE = 'plant', obj_, state))
 }

R/RcppR6.R

@@ -1,6 +1,6 @@
 ## Generated by RcppR6: do not edit by hand
 ## Version: 0.2.4
-## Hash: f6b511a10ca3bb03a599e71a6fa83092
+## Hash: 2149ebb5397c92102688df8a4c991aa7
 
 ##' @importFrom Rcpp evalCpp
 ##' @importFrom R6 R6Class
@@ -4814,6 +4814,30 @@ StochasticPatchRunner <- function(T, E) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
+      rainfall_init = function(x, y) {
+        FF16_Environment__rainfall_init(self, x, y)
+      },
+      rainfall_recompute = function() {
+        FF16_Environment__rainfall_recompute(self)
+      },
+      rainfall_eval = function(u) {
+        FF16_Environment__rainfall_eval(self, u)
+      },
+      rainfall_eval_range = function(u) {
+        FF16_Environment__rainfall_eval_range(self, u)
+      },
+      rainfall_add_point = function(x, y) {
+        FF16_Environment__rainfall_add_point(self, x, y)
+      },
+      rainfall_add_point_sorted = function(x, y) {
+        FF16_Environment__rainfall_add_point_sorted(self, x, y)
+      },
+      rainfall_add_points = function(x, y) {
+        FF16_Environment__rainfall_add_points(self, x, y)
+      },
+      rainfall_clear_points = function() {
+        FF16_Environment__rainfall_clear_points(self)
+      },
       canopy_openness = function(height) {
         FF16_Environment__canopy_openness(self, height)
       },
@@ -4823,9 +4847,6 @@ StochasticPatchRunner <- function(T, E) {
       set_fixed_environment = function(value, height_max) {
         FF16_Environment__set_fixed_environment(self, value, height_max)
       },
-      set_soil_infiltration_rate = function(rate) {
-        FF16_Environment__set_soil_infiltration_rate(self, rate)
-      },
       set_soil_water_state = function(state) {
         FF16_Environment__set_soil_water_state(self, state)
       },

R/ff16w.R

@@ -64,16 +64,17 @@ FF16w_StochasticPatchRunner <- function(p) {
 ## Helper:
 ##' @export
 ##' @rdname FF16_Environment
-##' @param infil_rate rate of water entering the first layer
-##' @param n_layers the number of layers
-##' @param init starting conditions
+##' @param soil_number_of_depths the number of soil layers
+##' @param rainfall_x list of x (time) positions for rainfall spline
+##' @param rainfall_y list of y (rainfall) positions for rainfall spline
 FF16w_make_environment <- function(canopy_light_tol = 1e-4, 
                                    canopy_light_nbase = 17,
                                    canopy_light_max_depth = 16, 
                                    canopy_rescale_usually = TRUE,
                                    soil_number_of_depths = 1,
                                    soil_initial_state = 0.0,
-                                   soil_infiltration_rate = 0.0) {
+                                   rainfall_x = seq(0, 9, 1),
+                                   rainfall_y = rep(0, 10)) {
 
   if(soil_number_of_depths < 1)
     stop("FF16w Environment must have at least one soil layer")
@@ -85,15 +86,15 @@ FF16w_make_environment <- function(canopy_light_tol = 1e-4,
                      canopy_light_nbase, 
                      canopy_light_max_depth)
 
-  # there might be a better way to skip this if using defaultss
+  # there might be a better way to skip this if using defaults
   if(sum(soil_initial_state) > 0.0) {
     if(soil_number_of_depths != length(soil_initial_state))
       stop("Not enough starting points for all layers")
 
     e$set_soil_water_state(soil_initial_state)
   }
 
-  e$set_soil_infiltration_rate(soil_infiltration_rate)
+  e$rainfall_init(rainfall_x, rainfall_y)
 
   return(e)
 }

inst/RcppR6_classes.yml

@@ -871,6 +871,30 @@ FF16_Environment:
       FF16_Environment object
       @export
   methods:
+    rainfall_init:
+      args: [x: "std::vector<double>", y: "std::vector<double>"]
+      return_type: void
+    rainfall_recompute:
+      args: []
+      return_type: void
+    rainfall_eval:
+      args: [u: double]
+      return_type: double
+    rainfall_eval_range:
+      args: [u: "std::vector<double>"]
+      return_type: "std::vector<double>"
+    rainfall_add_point:
+      args: [x: double, y: double]
+      return_type: void
+    rainfall_add_point_sorted:
+      args: [x: double, y: double]
+      return_type: void
+    rainfall_add_points:
+      args: [x: "std::vector<double>", y: "std::vector<double>"]
+      return_type: void
+    rainfall_clear_points:
+      args: []
+      return_type: void
     canopy_openness:
       args: [height: double]
       return_type: double
@@ -879,9 +903,6 @@ FF16_Environment:
     set_fixed_environment:
       args: [value: double, height_max: double]
       return_type: void
-    set_soil_infiltration_rate:
-      args: [rate: double]
-      return_type: void
     set_soil_water_state:
       args: [state: "std::vector<double>"]
       return_type: void

inst/include/plant/interpolator.h

@@ -16,16 +16,19 @@ class Interpolator {
   void initialise();
 
   void add_point(double xi, double yi);
+  void add_point_sorted(double xi, double yi);
   void clear();
 
   double eval(double u) const;
+  double operator()(double u) const;
   size_t size() const;
 
   double min() const;
   double max() const;
 
   std::vector<double> get_x() const;
   std::vector<double> get_y() const;
+
 
   // * R interface
   SEXP r_get_xy() const;

inst/include/plant/models/ff16_environment.h

@@ -4,6 +4,8 @@
 
 #include <plant/environment.h>
 #include <plant/canopy.h>
+#include <plant/interpolator.h>
+#include <iostream>
 
 using namespace Rcpp;
 
@@ -19,14 +21,59 @@ class FF16_Environment : public Environment {
       : canopy_rescale_usually(canopy_rescale_usually) {
     time = 0.0;
     canopy = Canopy();
-    soil_infiltration_rate = 0.0;
     vars = Internals(soil_number_of_depths);
     set_soil_water_state(std::vector<double>(soil_number_of_depths, 0.0));
+    rainfall_spline = interpolator::Interpolator();
   };
 
+  // first spline computation
+  void rainfall_init(std::vector<double> const& x, std::vector<double> const& y) {
+    rainfall_spline.init(x, y); // set x, y points and compute spline
+  }
+
+  // if any points are added to the spline, need to recompute spline
+  void rainfall_recompute() {
+    rainfall_spline.initialise();
+  }
+
+  double rainfall_eval(double u) const {
+    return rainfall_spline.eval(u);
+  }
+
+  std::vector<double> rainfall_eval_range(std::vector<double> const& u) const {
+    return rainfall_spline.r_eval(u);
+  }
+
+  /* if xi is not larger than every other x in the spline, recomputing the spline 
+  will throw an error and will need to clear_points() and reconstruct spline */
+  void rainfall_add_point(double xi, double yi) {
+    rainfall_spline.add_point(xi, yi);
+  }
+
+  // safe add point
+  void rainfall_add_point_sorted(double xi, double yi) {
+    rainfall_spline.add_point_sorted(xi, yi);
+  }
+
+  /* not sure if these two are necessary for the rainfall model. could also move these
+  to be general interpolator.h functions? */
+  void rainfall_add_points(std::vector<double> const& x, std::vector<double> const& y) {
+    util::check_length(x.size(), y.size());
+    for (size_t i = 0; i < x.size(); ++i) {
+      rainfall_add_point(x[i], y[i]); // refactor?
+    }
+  }
+
+  // rainfall_add_points_sorted ?
+
+  void rainfall_clear_points() {
+    rainfall_spline.clear();
+  }
+
   // Light interface
   bool canopy_rescale_usually;
 
+  // private?
   Canopy canopy;
 
   // Should this be here or in canopy?
@@ -51,12 +98,10 @@ class FF16_Environment : public Environment {
     canopy.r_init_interpolators(state);
   }
 
-  // Soil interface
-  double soil_infiltration_rate;
 
   virtual void compute_rates() {
     for (size_t i = 0; i < vars.state_size; i++) {
-      vars.set_rate(i, soil_infiltration_rate / (i+1));
+      vars.set_rate(i, rainfall_eval(time) / (i+1));
     }
   }
 
@@ -74,10 +119,6 @@ class FF16_Environment : public Environment {
     }
   }
 
-  void set_soil_infiltration_rate(double rate) {
-    soil_infiltration_rate = rate;
-  }
-
   // Core functions
   template <typename Function>
   void compute_environment(Function f_compute_competition, double height_max) {
@@ -92,7 +133,9 @@ class FF16_Environment : public Environment {
   void clear_environment() {
     canopy.clear();
   }
-
+
+private:
+  interpolator::Interpolator rainfall_spline;
 };
 
 inline Rcpp::NumericMatrix get_state(const FF16_Environment environment) {