Introduction

General flowchart
General flowchart

Workplan: main modules

access inventory data  
initialiazie stands list  
  
choose appropriate:  
    * Site Index classification tool  
    * height by dbh estimation model  
    * dbh growth function  
  
CHECK  
    all required info available?  
        YES, proceed to GROW  
    Age missing?  
        YES, goto (A)  
    Site Index missing?  
        NO, goto (A)  
        compute dominant height  
        assess Site Index  
        update stands list  
        (back to (0))  
    (A) other problems (under dev.), STOP  

GROW  
    (1) compile stand table  
    estimate dominant height and mean height for all dbh classes  
    [evaluate if thinning is required and execute treatment]   
    estimate dbh growth for each class  
    update stand table  
    ending condition attained?  
        NO, (back to (1))  
final stand table - produce summary

Implementation

access inventory data

library(googlesheets)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(tidyr)
library(magrittr)
## 
## Attaching package: 'magrittr'
## The following object is masked from 'package:tidyr':
## 
##     extract
getwd()
## [1] "C:/Users/ro/Documents/RICERCA/GitHub/DouglasFirGrowthModel_Italy/GO_DoNaTo"
source("Functions/access2DB.R", echo=F, verbose=F)
## Sheet successfully identified: "GO_baseDati"
# t <- 
tibble::tribble(
  ~table,              ~sheet,
  "species",           "Specie",
  "stands",            "AdS",
  "dbh_tally",         "Cavallettamento",
  "heights",           "Altezze",
  "increment_borings", "Carotine") %>% 
  group_by(sheet) %>% 
  do(assign(.$table, gs_read(DataBase, ws=.$sheet), pos=1)) %>%
  rm()
## Accessing worksheet titled 'AdS'.
## Parsed with column specification:
## cols(
##   id_ads = col_character(),
##   AdS_etichetta = col_character(),
##   quota = col_integer(),
##   `pend (%)` = col_integer(),
##   espos = col_character(),
##   eta2017 = col_integer(),
##   etaInizialeCarotina = col_integer(),
##   n.anelloPerEtàInizialeCarotina = col_integer(),
##   note = col_character(),
##   sup_m_quadri = col_integer(),
##   data_rilievo = col_character(),
##   NOTE_SELVICOLTURA = col_character(),
##   rilevatori = col_character()
## )
## Accessing worksheet titled 'Altezze'.
## Parsed with column specification:
## cols(
##   id_ads = col_character(),
##   sp = col_character(),
##   diam = col_integer(),
##   h = col_double(),
##   cld = col_integer()
## )
## Accessing worksheet titled 'Carotine'.
## Parsed with column specification:
## cols(
##   id_ads = col_character(),
##   id_fusto = col_integer(),
##   d130 = col_integer(),
##   `da-sottoc-a-5anni` = col_integer(),
##   `da-sottoc-a-CIRCA_40anniDiEtà` = col_integer(),
##   `(possibilmente senza mettere la virgola dei decimali)` = col_character(),
##   cld = col_integer()
## )
## Accessing worksheet titled 'Cavallettamento'.
## Parsed with column specification:
## cols(
##   id_ads = col_character(),
##   d_130 = col_integer(),
##   dou = col_integer(),
##   abe = col_integer(),
##   tig = col_integer(),
##   fag = col_integer(),
##   pin = col_integer(),
##   cas = col_integer(),
##   ont = col_integer(),
##   cer = col_integer()
## )
## Accessing worksheet titled 'Specie'.
## Parsed with column specification:
## cols(
##   sp = col_character(),
##   specie = col_character()
## )
# save(list=as.vector(t$table), file="GO_DonatoDB.Rdata")
# load("../GO_DonatoDB.Rdata")
dbh_tally <- dbh_tally %>% 
  gather("species", "freq", -id_ads, -d_130, na.rm=T) %>%
  left_join(select(stands, id_ads, sup_m_quadri)) %>%
  mutate(n_ha = 10000 * freq / sup_m_quadri) 
## Joining, by = "id_ads"
# dominant diameter
cum_frq <- 
  dbh_tally %>%
  group_by(id_ads, d_130) %>%
  summarise(n_ha_tot = sum(n_ha)) %>%
  group_by(id_ads) %>%
  arrange(desc(d_130)) %>%
  do(mutate(., cum_n_ha = cumsum(n_ha_tot)))
d_dom <-
  cum_frq %>%
  filter(cum_n_ha >= 100) %>%
  group_by(id_ads) %>%
  filter(d_130 == max(d_130)) %>%
  mutate(n_ha_tot = n_ha_tot - (cum_n_ha - 100)) %>%
  union(filter(cum_frq, cum_n_ha < 100)) %>%
  summarise(n_tot = sum(n_ha_tot), 
            d_dom = sqrt(sum(n_ha_tot * d_130^2)/n_tot ))

source("Functions/F_SiteIndex_MN94.R")
current_year <- 2017
# dominant height and 
# site class  evaluation
h_dom <-
  heights %>%
  filter(sp == "dou") %>%
  group_by(id_ads) %>%
  do(hm=loess(h~diam, .)) %>%
  left_join(d_dom) %>%
  group_by(id_ads) %>%
  do(d_dom = .$d_dom, 
     f = predict(.$hm[[1]], data.frame(diam=.$d_dom), se = TRUE)) %>%
  group_by(id_ads) %>%
  do(d_dom = .$d_dom[[1]], 
     h_dom = (.$f[[1]]$fit[[1]]),
     se.h_dom = (.$f[[1]]$se.fit[[1]])
     ) %>%
  left_join(select(stands, id_ads, eta2017)) %>%
  mutate(plantation_year = 2017 - eta2017) %>%
  select(-eta2017) %>%
  mutate(SI_MN94_50 = dougSI_MN94(current_year - plantation_year, h_dom))
## Joining, by = "id_ads"
## Warning: Grouping rowwise data frame strips rowwise nature

## Warning: Grouping rowwise data frame strips rowwise nature
## Joining, by = "id_ads"

[[to be continued …]]