More MSR/aspatial fixes

Hopefully MSR should be working correctly now, both with single value and with map as rules input

Author: wburke24

R/CreateFlownet.R

@@ -5,17 +5,18 @@
 #' @param readin readin indicates the maps to be used. If CreateFlowmet.R is run it's own, this should point to the template which references the maps(or values)
 #' used for the various levels and statevars. Otherwise, if run inside of RHESSysPreprocess, readin will use the map data from world_gen.R,
 #' Streams map, and other optional maps, still need to be specified.
-#' @param asp_rules List of aspatial structure and inputs. Also can be path to rules file - must be used along with template input.
+#' @param asprules List of aspatial structure and inputs. Also can be path to rules file - must be used along with template input.
 #' @param road_width >0, defaults to 1.
 #' @inheritParams RHESSysPreprocess
 #' @author Will Burke
 #' @export
 
 CreateFlownet = function(flownet_name,
                          readin = NULL,
+                         template = NULL,
                          type = "raster",
                          typepars = NULL,
-                         asp_rules = NULL,
+                         asprules = NULL,
                          streams = NULL,
                          overwrite = FALSE,
                          roads = NULL,
@@ -39,11 +40,12 @@ CreateFlownet = function(flownet_name,
   if (!is.logical(overwrite)) {stop("overwrite must be logical")} # check overwrite inputs
   if (file.exists(flownet_name) & overwrite == FALSE) {stop(noquote(paste("Flowtable",flownet_name,"already exists.")))}
 
+  template_list = template_read(template)
+  template_clean = template_list[[1]] # template in list form
+  var_names = template_list[[2]] # names of template vars
+  map_info = template_list[[5]] # tables of maps and their inputs/names in the template
+
   if (!wrapper & is.character(readin)) { #if run outside of rhessyspreprocess.R, and if readin is character. readin is the template (and path)
-    template_list = template_read(readin)
-    template_clean = template_list[[1]] # template in list form
-    var_names = template_list[[2]] # names of template vars
-    map_info = template_list[[5]] # tables of maps and their inputs/names in the template
     cfmaps = rbind(map_info,c("cell_length","none"), c("streams","none"), c("roads","none"), c("impervious","none"),c("roofs","none"))
   } else if (wrapper | (!wrapper & is.matrix(readin))) { # map info is passsed directly from world gen - either in wrapper or outside of wrapper and readin is matrix
     cfmaps = readin
@@ -116,22 +118,28 @@ CreateFlownet = function(flownet_name,
   raw_impervious_data = NULL
   if (!is.null(impervious)) {raw_impervious_data =  map_list[[cfmaps[cfmaps[,1] == "impervious",2]]]}
 
+  # TODO - this should get cleaned up eventually - there's at least 2 different versions of the maps being used
+  if (length(readmap@data[,1]) == 1) {
+    map_df = as.data.frame(readmap@data) # works for 1 patch world
+  } else {
+    map_df = as.data.frame(readmap) #make data frame for ease of use
+  }
+
   # read aspatial rules if needed
-  if (!is.null(asp_rules) && is.character(asp_rules)) {
+  if (!is.null(asprules)) {
     asp_map = template_clean[[which(var_names == "asp_rule")]][3] # get rule map/value
-    if (suppressWarnings(any(is.na(as.numeric(asp_map))))) { # if it's a map
+
+    if (!"asp_rule" %in% notamap) { # if it's a map
       asp_map = gsub(".tif|.tiff","",asp_map)
       asp_mapdata = as.data.frame(readmap)[asp_map]
 
       # --- doing manipulation of the asp map if needed here ---
-      if (template_clean[[which(var_names == "asp_rule")]][2] == "mode" ) {
-
+      if (template_clean[[which(var_names == "asp_rule")]][2] == "mode") {
         # mode for aggregating by mode
         mode_fun = function(x) {
           ux <- unique(x)
           ux[which.max(tabulate(match(x, ux)))]
         }
-
         # fun to check if rules are all unique per patch (kinda slow but oh well) -----
         check_rules = function(patches, asp_rules) {
           patches_u = unique(patches[!is.na(patches)])
@@ -155,28 +163,22 @@ CreateFlownet = function(flownet_name,
                                map_info[map_info[, 1] == "strata", 2]))
         level_names = unique(gsub(".tiff|.tif|.asc","",level_names))
 
-        if (length(readmap@data[,1]) == 1) {
-          map_df = as.data.frame(readmap@data) # works for 1 patch world
-        } else {
-          map_df = as.data.frame(readmap) #make data frame for ease of use
-        }
-
-        asp_maps = aggregate(map_df$asprule, by = map_df[level_names], FUN = mode_fun, simplify = T)
+        asp_maps = aggregate(map_df[[cfmaps[cfmaps[,1] == "asp_rule",2]]], by = map_df[level_names], FUN = mode_fun, simplify = T)
         names(asp_maps)[which(names(asp_maps) == "x")] = "asprule"
 
         rules_out = check_rules(patches = asp_maps$pch_30m1000, asp_rules = asp_maps$asprule)
         if (!is.null(rules_out)) {
           print(rules_out)
           stop("Mode of rules was attempted but there are still patches with multiple rules. Check input maps.")
         }
-
         asp_mapdata = asp_maps$asprule
       }
 
-    } else if (suppressWarnings(all(!is.na(as.numeric(asp_map))))) { # if is a single number
-      asp_mapdata = as.numeric(asp_map)
+    } else if ("asp_rule" %in% notamap) { # if is a single number
+      asp_mapdata = map_df[[cfmaps[cfmaps[,1] == "basin",2]]]
+      asp_mapdata[!is.na(asp_mapdata)] = as.numeric(cfmaps[cfmaps[,1] == "asp_rule",2])
     }
-    asp_list = aspatial_patches(asprules = asp_rules, asp_mapdata = asp_mapdata)
+    asp_list = aspatial_patches(asprules = asprules, asp_mapdata = asp_mapdata)
   }
 
   # ------------------------------ Make flownet list ------------------------------
@@ -197,15 +199,9 @@ CreateFlownet = function(flownet_name,
     skip_hillslope_check = skip_hillslope_check)
 
   # ------------------------------ Multiscale routing/aspatial patches ------------------------------
-  if (!is.null(asp_rules)) {
-    if ("asp_rule" %in% notamap) {
-      map_list[["asp_rule"]] = raw_basin_data
-      map_list[["asp_rule"]][!is.na(map_list[["asp_rule"]])] = as.numeric(cfmaps[cfmaps[,1] == "asp_rule",2])
-      #cfmaps = rbind(cfmaps, c("asp_rule", "asp_rule"))
-      cfmaps[cfmaps[,1] == "asp_rule", 2] = "asp_rule"
-    }
-
-    CF1 = multiscale_flow(CF1 = CF1, asp_maps = asp_maps, cfmaps = cfmaps, asp_list = asp_list)
+  if (!is.null(asprules)) {
+    patch_map = map_df[[cfmaps[cfmaps[,1] == "patch",2]]]
+    CF1 = multiscale_flow(CF1 = CF1, asp_map = asp_mapdata, patch_map = patch_map, asp_list = asp_list)
   }
 
   # ---------- Flownet list to flow table file ----------

R/RHESSysPreprocess.R

@@ -127,7 +127,7 @@ RHESSysPreprocess = function(template,
                             asprules = asprules)
 
   readin = world_gen_out[[1]]
-  asp_rules = world_gen_out[[2]]
+  #asp_rules = world_gen_out[[2]]
 
   # ---------- Run CreateFlownet ----------
   cat("Begin CreateFlownet.R")
@@ -143,9 +143,10 @@ RHESSysPreprocess = function(template,
 
   CreateFlownet(flownet_name = flownet_name,
                 readin = readin,
+                template = template,
                 type = type,
                 typepars = typepars,
-                asp_rules = asp_rules,
+                asprules = asprules,
                 streams = streams,
                 overwrite = overwrite,
                 roads = roads,

R/aspatial_patches.R

@@ -3,7 +3,7 @@
 #' @param asp_mapdata map data or value, indicating the rule IDs being used
 #' @author Will Burke
 
-aspatial_patches = function(asprules,asp_mapdata) {
+aspatial_patches = function(asprules, asp_mapdata) {
 
   # some functions
   # splits a vector up at a specific character
@@ -49,7 +49,8 @@ aspatial_patches = function(asprules,asp_mapdata) {
 
   # ---------- build rulevars based on rules and map data ----------
   map_ids = unique(asp_mapdata) # get rule IDs from map/input
-  map_id_tags = paste("rule_",map_ids,sep = "") # all map IDs concated w tags for referencing/reading in code
+  if (is.list(map_ids)) {map_ids = unlist(map_ids)}
+  map_id_tags = paste0("rule_", as.numeric(map_ids)) # all map IDs concated w tags for referencing/reading in code
 
   asp_vars = as.list(rep(0,length(map_ids))) # highest level list of the different rules
   # strata_index = as.list(rep(0,length(map_ids))) # get rid of in this version i think

R/multiscale_flow.R

@@ -2,12 +2,12 @@
 #'
 #' Modifies an existing flowtable (R list) to work with multisacle routing.
 #' @param CF1 Flowtable list object created from make_flow_list
-#' @param map_list List containing input maps as matrices
-#' @param cfmaps Table of map types and file names
+#' @param asp_map Map of the aspatial rules
+#' @param patch_map Map of the patches
 #' @param asp_list List of aspatial rules
 # Will Burke 1/16/19
 
-multiscale_flow = function(CF1, asp_maps, cfmaps, asp_list) {
+multiscale_flow = function(CF1, asp_map, patch_map, asp_list) {
 
   #nbr 121
   #nbr patch 11955
@@ -53,11 +53,9 @@ multiscale_flow = function(CF1, asp_maps, cfmaps, asp_list) {
   cat("Creating multiscale flowtable - this may take a moment with many patches")
 
   # ----- Variable setup -----
-  asp_map = asp_maps[[cfmaps[cfmaps[, 1] == "asp_rule", 2]]] # matrix of aspatial rules
-
   patch_ID = unlist(lapply(CF1, "[[", 9)) # patch IDs from cf1
   numbers = unlist(lapply(CF1, "[[", 1)) # flow list numbers
-  raw_patch_data = asp_maps[[cfmaps[cfmaps[, 1] == "patch", 2]]] # get patch matrix inside the function
+  raw_patch_data = patch_map
   rulevars = asp_list # get rules - state variable overrides
   CF2 = list() # empty list for new flow list
 

R/world_gen.R

@@ -114,7 +114,7 @@ world_gen = function(template,
     rulevars = aspatial_patches(asprules = asprules, asp_mapdata = asp_mapdata)
 
     if (is.data.frame(asp_mapdata)) { # add ruleID to levels df
-      asp_level = aggregate(asp_mapdata$asprule.tiff, by = as.list(levels), FUN = mode_fun)
+      asp_level = aggregate(asp_mapdata[[1]] , by = as.list(levels), FUN = mode_fun)
       levels = merge(levels, asp_level[,c("h", "z", "p", "x")], by = c("h", "z", "p"), sort = F)
       levels = levels[c("w", "b", "h", "z", "p", "s", "x")]
       names(levels)[7] = "a"