additional fixes for aspatial patches

Author: wburke24

R/CreateFlownet.R

@@ -122,10 +122,61 @@ CreateFlownet = function(flownet_name,
     if (suppressWarnings(any(is.na(as.numeric(asp_map))))) { # 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" ) {
+
+        # 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)])
+          patches_i = lapply(patches_u, function(x, y) {which(x == y)}, patches)
+          outrules = lapply(patches_i, function(x, y) {y[x]}, asp_rules)
+          lens = sapply(outrules, function(x) {length(unique(x))})
+          lens_i = which(unlist(lens) > 1)
+          if (length(lens_i) > 0) {
+            cat("There are patches with multiple multiscale rules:\n")
+            out = data.frame(patch = patches[lens_i], rules = sapply(outrules[lens_i], function(x) {paste(unique(x), collapse=", ") }) )
+            return(out)
+          } else {
+            return(NULL)
+          }
+        }
+
+        level_names = unname(c(map_info[map_info[, 1] == "basin", 2],
+                               map_info[map_info[, 1] == "hillslope", 2],
+                               map_info[map_info[, 1] == "zone", 2],
+                               map_info[map_info[, 1] == "patch", 2],
+                               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)
+        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)
     }
-    asp_rules = aspatial_patches(asprules = asp_rules, asp_mapdata = asp_mapdata)
+    asp_list = aspatial_patches(asprules = asp_rules, asp_mapdata = asp_mapdata)
   }
 
   # ------------------------------ Make flownet list ------------------------------
@@ -154,7 +205,7 @@ CreateFlownet = function(flownet_name,
       cfmaps[cfmaps[,1] == "asp_rule", 2] = "asp_rule"
     }
 
-    CF1 = multiscale_flow(CF1 = CF1, map_list = map_list, cfmaps = cfmaps, asp_list = asp_rules)
+    CF1 = multiscale_flow(CF1 = CF1, asp_maps = asp_maps, cfmaps = cfmaps, asp_list = asp_list)
   }
 
   # ---------- Flownet list to flow table file ----------

R/RHESSysPreprocess.R

@@ -153,7 +153,8 @@ RHESSysPreprocess = function(template,
                 roofs = roofs,
                 wrapper = wrapper,
                 parallel = parallel,
-                make_stream = make_stream)
+                make_stream = make_stream,
+                skip_hillslope_check = skip_hillslope_check)
 
   # ---------- Run build_meta ----------
   # if (meta) {

R/aspatial_patches.R

@@ -48,7 +48,7 @@ aspatial_patches = function(asprules,asp_mapdata) {
   # strata data frames (num = patch ct) - rows are statevars, cols are strata
 
   # ---------- build rulevars based on rules and map data ----------
-  map_ids = unique(asp_mapdata)[[1]] # get rule IDs from map/input
+  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
 
   asp_vars = as.list(rep(0,length(map_ids))) # highest level list of the different rules

R/multiscale_flow.R

@@ -7,18 +7,23 @@
 #' @param asp_list List of aspatial rules
 # Will Burke 1/16/19
 
-multiscale_flow = function(CF1, map_list, cfmaps, asp_list) {
+multiscale_flow = function(CF1, asp_maps, cfmaps, asp_list) {
+
+  #nbr 121
+  #nbr patch 11955
+  # no times nbr asp ct is empty
+  #which(lapply(CF1,"[[", 9) == 15660)
 
   # functions for applys
   apply_patches = function(CFp) {
     id = paste0("rule_", unique(asp_map[which(raw_patch_data == CFp$PatchID)]))
     asp_count = ncol(rulevars[[id]]$patch_level_vars[1, ]) - 1 # get number of aspatial patches for current patch
     asp = c(1:asp_count)
-    CFasp = lapply(asp, add_asp,CFp)
+    CFasp = lapply(asp, add_asp,CFp, id)
     unlist(CFasp,recursive = F)
     return(CFasp)
   }
-  add_asp = function(asp, CFp) {
+  add_asp = function(asp, CFp, id) {
     CFasp = CFp
     CFasp$PatchID = CFp$PatchID * 100 + asp # aspatial patch ID is old patch ID *100 + aspatial number
     CFasp$Number = CFp$Number * 100 + asp # same modification to number
@@ -48,11 +53,11 @@ multiscale_flow = function(CF1, map_list, cfmaps, asp_list) {
   cat("Creating multiscale flowtable - this may take a moment with many patches")
 
   # ----- Variable setup -----
-  asp_map = map_list[[cfmaps[cfmaps[, 1] == "asp_rule", 2]]] # matrix of aspatial rules
+  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 = map_list[[cfmaps[cfmaps[, 1] == "patch", 2]]] # get patch matrix inside the function
+  raw_patch_data = asp_maps[[cfmaps[cfmaps[, 1] == "patch", 2]]] # get patch matrix inside the function
   rulevars = asp_list # get rules - state variable overrides
   CF2 = list() # empty list for new flow list
 

R/world_gen.R

@@ -96,6 +96,11 @@ world_gen = function(template,
   levels = as.data.frame(levels)
   colnames(levels) = c("w", "b", "h", "z", "p", "s")
 
+  # mode for aggregating by mode
+  mode_fun = function(x) {
+    ux <- unique(x)
+    ux[which.max(tabulate(match(x, ux)))]
+  }
 
   # -------------------- Aspatial Patch Processing --------------------
   rulevars = NULL
@@ -108,8 +113,12 @@ world_gen = function(template,
     }
     rulevars = aspatial_patches(asprules = asprules, asp_mapdata = asp_mapdata)
 
-    if (is.data.frame(asp_mapdata)) { # add ruleID to levels matrix
-      levels = cbind(levels, a = unname(as.matrix(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)
+      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"
+      #levels = cbind(levels, a = unname(as.matrix(asp_mapdata)))
     } else if (is.numeric(asp_mapdata)) {
       levels = cbind(levels, a = rep(asp_mapdata,length(levels[,1])) )
     }
@@ -154,11 +163,7 @@ world_gen = function(template,
   } # end loop through n_basestations
 
   # -------------------- Process Template + Maps --------------------
-  # mode for aggregating by mode
-  mode_fun = function(x) {
-    ux <- unique(x)
-    ux[which.max(tabulate(match(x, ux)))]
-  }
+
   # Build list based on operations called for by template
   statevars = vector("list",length(template_clean))