R/joinCWDData.R

In KateMMiller/forestMIDNarch: Package for importing, joining and querying MS Access based MIDN forest database tables.

#' @include joinLocEvent.R
#'
#' @importFrom dplyr select filter arrange mutate summarise group_by case_when
#' @importFrom magrittr %>%
#'
#' @title joinCWDData: compile coarse woody debris data.
#'
#' @description This function combines calculates CWD volume for each plot. Must run importData first.
#' Note that plots in Deer Exclosures are not sampled for CWD and will not be returned.
#'
#' @param units Calculates CWD Volume based on different units.
#' \describe{
#' \item{"ha"}{Default. Returns CWD volume as cubic m/hectare}
#' \item{"acres"}{Returns CWD volume as cubic ft/ acre}
#'}
#'
#' @return returns a dataframe with CWD volume for each plot, one with cubic m/ha and cubic ft/acre
#'
#' @examples
#' importData() #imports using default odbc
#' # Compile CWD data for FRSP for most recent survey and return in ft^3/acre
#' cwd_data <- joinCWDData(park = 'FRSP', from = 2015, to = 2018, units = 'acres')
#'
#' # Compile CWD data for all parks (default) for most recent survey and return in m^3/ha (default)
#' cwd_data <- joinCWDData(from = 2015, to = 2018)
#'
#' @export
#'
#------------------------
# Join CWD table and filters by park, year, and plot/visit type
#------------------------
joinCWDData<-function(units=c('ha','acres'), park='all',from=2007, to=2019, QAQC=FALSE, locType='VS', output, ...){
  units<-match.arg(units)
  # Prepare the CWD data
  park.plots <- force(joinLocEvent(park=park, from=from,to=to,QAQC=QAQC,locType=locType, output='verbose'))
  park.plots <- park.plots %>% filter(Loc_Type!='Deer') %>%
    select(Location_ID, Event_ID, Unit_Code, Plot_Name, Plot_Number, X_Coord, Y_Coord, Panel, Year, Event_QAQC,
           cycle) %>% droplevels()

  cwd1<-merge(park.plots,cwd,by='Event_ID', all.x=T,all.y=F)
  cwd2<-merge(cwd1[,c("Event_ID","TSN","Diameter","Decay_Class_ID","Hollow","Transect","Distance","Wood_Type")],
              plants[,c('TSN','Latin_Name','Common')], by='TSN',all.x=T)
  cwd.std<-merge(cwd2,stand[,c("Event_ID","Slope_UP","Slope_BR","Slope_BL")],by="Event_ID",all.x=T,all.y=F)

  cwd.std2<-cwd.std %>% mutate(slope=case_when(Transect=='BL'~ Slope_BL, Transect=='BR'~ Slope_BR, Transect=='UP'~ Slope_UP),
                               pct.slope=ifelse(is.na(slope),0,tan(slope*pi/180)*100),
                               hdist=((((pct.slope/100)^2)+1)^0.5)*((pi^2)/(8*15)),
                               diam=Diameter^2)
  cwd3<-cwd.std2 %>% group_by(Event_ID,Transect, hdist, Latin_Name, Decay_Class_ID) %>% summarise(diam=sum(diam)) # Add species and decay class
  cwd4<-cwd3 %>% group_by(Event_ID, Latin_Name, Decay_Class_ID) %>% summarise(CWD_Vol=ifelse(is.na(sum(diam)),0,sum(hdist*diam)/3))
  cwd5<-merge(park.plots,cwd4,by="Event_ID", all.x=T)

  cwd6<-if (units=='acres'){
    cwd5 %>% mutate(CWD_Vol=CWD_Vol*35.314667/2.4710538)
    # 35.314667 is the # cubic feet in a cubic meter. 2.4710538 is # acres in 1 hectare.)
  } else if (units=='ha'){return(cwd5)
  }

  cwd7<-merge(park.plots,cwd6[,c("Event_ID","CWD_Vol", "Latin_Name","Decay_Class_ID")],by="Event_ID",all.x=T)
  cwd7[,"CWD_Vol"][is.na(cwd7[,"CWD_Vol"])]<-0

  return(data.frame(cwd7))
} # end of function