flukeclimatemse: Run Climate Fluke MSE

# This script contains the operating model function: dosurveys

# ToDo:
  # Search for #Letter equation numbers & data table sources
  # Update survey catchability using bottomtrawl survey from NOAA

#' @title Generate an index of biomass at age
#'
#' @description This function uses stock abundance at age to generate an index of biomass at age for use as the stock assessment survey data, assuming that catchability at age for the survey remains constant over time.
#'
#' @param N_abund A vector of abundance at age at the start of the year (t), no default.
#'
#' @return A list containing:
#' Survey biomass in year t and corresponding survey observation error.
#'
#' @family operating model functions
#'
#' @examples
#' # Example
#' dosurveys(N_abund = c(12000, 13000, 12900, 12800, 12500, 12000, 11000, 10000)) # fish abundance for age 0-7+
#'

dosurveys <- function(N_abund = NULL){
  # Fixed/locally calculated inputs
  q_age <- c(0.06, 0.06, 0.42, 0.66, 0.70, 0.70, 0.59, 0.58) # In Table #H ??? Catchability at age 0-7+ for survey index, see data_calculations.R for details, pulled from 66th Northeast Regional Stock Assessment Workshop Assessment Report (NEFSC 2019)

  # Calculate survey observation error
  surveyError <- rnorm(1,mean = 0, sd = 0.178) # Equation #K ???, standard deviation set to average of historic spring & fall survey CV see data_calculations.R for details
  # Calculate survey index (index of biomass at age) in year t
  survey <- sum(q_age*N_abund) + surveyError # Equation #J ???

  return(list(OM_Survey=survey, # Survey biomass in year t
              OM_surveyError = surveyError)) # Observation error of survey biomass in year t
}




# # Uncomment to run examples
#
#
# # Example abundance data
# # 2018 assessment abundance estimates at age from 66th Northeast Regional Stock Assessment Workshop Assessment Report (NEFSC 2019) Table A89, pg240 2018, years 1982-2017
# AbundData <- read.csv("/Users/ahart2/Research/Summer_Flounder_MSE/Data/SourceFiles/Abund_at_age_A89_NEFSC66AssessmentReport.csv", header=TRUE)
# AbundData <- AbundData[,2:9]
# for(icol in 1:ncol(AbundData)){
#   AbundData[,icol] <- as.numeric(gsub(",","",AbundData[,icol]))
# }
# colnames(AbundData) <- c("Age0", "Age1", "Age2", "Age3", "Age4", "Age5", "Age6", "Age7plus")
#
# # Example
# dosurveys(N_abund = unlist(AbundData[nrow(AbundData),]))