Commit de9a9a95 authored by bberges's avatar bberges
Browse files

update repo with new scripts after post-cruise meeting

parent ea790079
This diff is collapsed.
# script that derives the HERAS index directly from the StoX projects.
rm(list=ls())
library(Rstox)
library(FLFleet)
library(ggplotFL)
library(mgcv)
library(tidyverse)
library(ggplot2)
#Set up directories
#path <- 'J:/git/HERAS/'
path <- 'J:/git/heras_index_kbwot'
try(setwd(path),silent=TRUE)
mainPath <- file.path(".")
dataPath <- file.path(".","data")
stoxPath <- file.path(".","StoX")
outPath <- file.path(".","output")
functionPath <- file.path(".","functions")
resultsPath <- file.path(".","results")
figurePath <- file.path(".","figures")
load(file.path(resultsPath,'HERAS_baseline.RData'))
load(file.path(resultsPath,'HERAS_sensitivity_split_genetics.RData'))
######################################
# compare NSAS
######################################
print(as.numeric(HERAS_sensitivity_split_genetics$stk.NSAS@stock[,ac(2020)]))
print(as.numeric(HERAS_baseline$stk.NSAS@stock[,ac(2020)]))
print(as.numeric((HERAS_baseline$stk.NSAS@stock[,ac(2020)]-HERAS_sensitivity_split_genetics$stk.NSAS@stock[,ac(2020)])/HERAS_baseline$stk.NSAS@stock[,ac(2020)]))
df.HERAS_genetics_strata.NSAS <- as.data.frame(HERAS_sensitivity_split_genetics$HERAS.NSAS$HERAS_strata@index[,ac(2020)])
df.HERAS_baseline_strata.NSAS <- as.data.frame(HERAS_baseline$HERAS.NSAS$HERAS_strata@index[,ac(2020)])
df.HERAS_genetics_all.NSAS <- as.data.frame(HERAS_sensitivity_split_genetics$HERAS.NSAS$HERAS_all@index[,ac(2020)])
df.HERAS_baseline_all.NSAS <- as.data.frame(HERAS_baseline$HERAS.NSAS$HERAS_all@index[,ac(2020)])
png(file.path(figurePath,'abu_baseline_genetics_strata_all_NSAS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_all.NSAS[df.HERAS_baseline_all.NSAS$season=='all',],aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_all.NSAS[df.HERAS_baseline_all.NSAS$season=='all',],aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age NSAS - all strata')
dev.off()
png(file.path(figurePath,'abu_baseline_genetics_strata_141_NSAS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_strata.NSAS[df.HERAS_baseline_strata.NSAS$season=='all' & df.HERAS_baseline_strata.NSAS$area == 141,],
aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_strata.NSAS[df.HERAS_genetics_strata.NSAS$season=='all' & df.HERAS_genetics_strata.NSAS$area == 141,],
aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age NSAS - strata 141')
dev.off()
png(file.path(figurePath,'abu_baseline_genetics_strata_11_NSAS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_strata.NSAS[df.HERAS_baseline_strata.NSAS$season=='all' & df.HERAS_baseline_strata.NSAS$area == 11,],
aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_strata.NSAS[df.HERAS_genetics_strata.NSAS$season=='all' & df.HERAS_genetics_strata.NSAS$area == 11,],
aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age NSAS - strata 11')
dev.off()
######################################
# compare WBSS
######################################
print(as.numeric(HERAS_sensitivity_split_genetics$stk.WBSS@stock[,ac(2020)]))
print(as.numeric(HERAS_baseline$stk.WBSS@stock[,ac(2020)]))
print(as.numeric((HERAS_baseline$stk.WBSS@stock[,ac(2020)]-HERAS_sensitivity_split_genetics$stk.WBSS@stock[,ac(2020)])/HERAS_baseline$stk.WBSS@stock[,ac(2020)]))
df.HERAS_genetics_strata.WBSS <- as.data.frame(HERAS_sensitivity_split_genetics$HERAS.WBSS$HERAS_strata@index[,ac(2020)])
df.HERAS_baseline_strata.WBSS <- as.data.frame(HERAS_baseline$HERAS.WBSS$HERAS_strata@index[,ac(2020)])
df.HERAS_genetics_all.WBSS <- as.data.frame(HERAS_sensitivity_split_genetics$HERAS.WBSS$HERAS_all@index[,ac(2020)])
df.HERAS_baseline_all.WBSS <- as.data.frame(HERAS_baseline$HERAS.WBSS$HERAS_all@index[,ac(2020)])
png(file.path(figurePath,'abu_baseline_genetics_strata_all_WBSS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_all.WBSS[df.HERAS_baseline_all.WBSS$season=='all',],aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_all.WBSS[df.HERAS_baseline_all.WBSS$season=='all',],aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age WBSS - all strata')
dev.off()
png(file.path(figurePath,'abu_baseline_genetics_strata_141_WBSS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_strata.WBSS[df.HERAS_baseline_strata.WBSS$season=='all' & df.HERAS_baseline_strata.WBSS$area == 141,],
aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_strata.WBSS[df.HERAS_genetics_strata.WBSS$season=='all' & df.HERAS_genetics_strata.WBSS$area == 141,],
aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age WBSS - strata 141')
dev.off()
png(file.path(figurePath,'abu_baseline_genetics_strata_11_WBSS.png'), width = 10, height = 10, units = "cm", res = 300, pointsize = 10)
ggplot()+
geom_point(data = df.HERAS_baseline_strata.WBSS[df.HERAS_baseline_strata.WBSS$season=='all' & df.HERAS_baseline_strata.WBSS$area == 11,],
aes(x=as.factor(age),y=data),col='red')+
geom_point(data = df.HERAS_genetics_strata.WBSS[df.HERAS_genetics_strata.WBSS$season=='all' & df.HERAS_genetics_strata.WBSS$area == 11,],
aes(x=as.factor(age),y=data),col='blue')+
xlab('age')+
ylab('Abundance')+
ggtitle('Abundance at age WBSS - strata 11')
dev.off()
###############################################
# dump
###############################################
load(file.path(resultsPath,'HERAS_strata_definition.RData'))
newStrata_HERAS.NSAS <- HERAS.NSAS
newStrata_HERAS.WBSS <- HERAS.WBSS
load(file.path(resultsPath,'HERAS_sensitivity_SA_error.RData'))
SA_error_HERAS.NSAS <- HERAS.NSAS
SA_error_HERAS.WBSS <- HERAS.WBSS
# plotting
#df_EU <- as.data.frame(baseline_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index)
#df_EU$run <- 'baseline'
newStrata_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index <- (newStrata_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index-baseline_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index)/baseline_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index
newStrata_df_EU <- as.data.frame(newStrata_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index)
newStrata_df_EU$run <- 'newStrata'
SA_error_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index <- (SA_error_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index-baseline_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index)/baseline_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index
SA_error_df_EU <- as.data.frame(SA_error_HERAS.NSAS$HERAS_EU_NO[,ac(2017:2019),,,'EU']@index)
SA_error_df_EU$run <- 'SA_error'
#df_all <- rbind(df_EU,newStrata_df_EU,SA_error_df_EU)
df_all <- rbind(newStrata_df_EU,SA_error_df_EU)
df_all$year <- factor(df_all$year)
df_all$age <- factor(df_all$age)
ggplot(df_all,aes(x=age,y=data,col=year))+
geom_point()+
facet_wrap(~run,scales = "free_y")
......@@ -10,6 +10,7 @@ library(gdata) # drop.levels
library(reshape) # melt and cast
library(gstat)
library(sp)
library(tidyr)
path <- 'J:/git/heras_index_kbwot'
......@@ -34,7 +35,7 @@ mkPlot <- TRUE
# load species list - in git data path
speciesList <- read.csv(file.path(dataPath,'species_codes_201911.csv'), fill = TRUE, header = TRUE)
surveyYearMat <- c(2019,2018,2017,2016,2015)
surveyYearMat <- c(2020,2019,2018,2017,2016,2015)
for(surveyYear in surveyYearMat){
reportPath <- file.path(".","reports",surveyYear)
......@@ -77,7 +78,6 @@ for(surveyYear in surveyYearMat){
dir.create(file.path(figurePath,Cruise$CruiseCountry),showWarnings = FALSE)
if(mkPlot){
source(file.path(functionPath,"plot_biotic.R"))
### plot Biotic data
plot_biotic(Cruise,
Biology,
......@@ -121,6 +121,62 @@ for(surveyYear in surveyYearMat){
}
}
BiologySelect <- BiologyAll[ BiologyAll$BiologyIndividualAge != '' &
BiologyAll$BiologyLengthClass != '' &
BiologyAll$CatchSpeciesCode == 126425,]
BiologySelect$BiologyIndividualAge <- as.numeric(as.character(BiologySelect$BiologyIndividualAge))
BiologySelect$BiologyLengthClass <- as.numeric(as.character(BiologySelect$BiologyLengthClass))
allData <- as.data.frame(cbind(BiologySelect$BiologyLengthClass,BiologySelect$BiologyIndividualAge))
colnames(allData) <- c('length','age')
ggplot()+
geom_point(data=allData,
aes(x=length,y=age),col='grey',alpha=0.1)+
geom_point(data=BiologySelect,
aes(x=BiologyLengthClass,y=BiologyIndividualAge),col='red')+
facet_wrap(~CruiseLocalID)
catchNL <- CatchAll[CatchAll$CruiseLocalID == 'NLHERAS2020',]
catchNL$CatchSpeciesCategoryWeight <- as.numeric(as.character(catchNL$CatchSpeciesCategoryWeight))
catchNL <- catchNL[catchNL$CatchSpeciesCategoryWeight > 1,]
catchNL$HaulNumber <- as.numeric(as.character(catchNL$HaulNumber))
catchNL$CatchSpeciesCode <- as.numeric(as.character(catchNL$CatchSpeciesCode))
catchNL$CatchSpeciesCategoryWeight <- as.numeric(as.character(catchNL$CatchSpeciesCategoryWeight))
catchNL <- as.data.frame(cbind(catchNL$HaulNumber,catchNL$CatchSpeciesCode,catchNL$CatchSpeciesCategoryWeight))
colnames(catchNL) <- c('HaulNumber','species','catch')
catchNL <- dcast(catchNL, HaulNumber~species,value.var = "catch",sum)
speciesNames <- as.character(speciesList$SPECIESNAME[match(as.numeric(colnames(catchNL[2:dim(catchNL)[2]])),speciesList$WORMS)])
colnames(catchNL)[2:dim(catchNL)[2]] <- speciesNames
catchNL$total <- rowSums(catchNL[,2:dim(catchNL)[2]])
#catchNL[,2:(dim(catchNL)[2]-1)] <- catchNL[,2:(dim(catchNL)[2]-1)]/apply(catchNL[,2:(dim(catchNL)[2]-1)],2,max)
HaulNL <- HaulAll[HaulAll$CruiseLocalID == 'NLHERAS2020',]
HaulNL$HaulNumber <- as.numeric(as.character(HaulNL$HaulNumber))
HaulNL$HaulStartLongitude <- as.numeric(as.character(HaulNL$HaulStartLongitude))
HaulNL$HaulStartLatitude <- as.numeric(as.character(HaulNL$HaulStartLatitude))
HaulNL <- as.data.frame(cbind(HaulNL$HaulNumber,HaulNL$HaulStartLatitude,HaulNL$HaulStartLongitude))
colnames(HaulNL) <- c('HaulNumber','lat','lon')
catchNL <- catchNL %>%
left_join(HaulNL)
ggplot(world, aes(long, lat)) +
geom_map(map=world, aes(map_id=region), fill=NA, color="black") +
coord_quickmap(xlim=c(floor(min(catchNL$lon)-0.5),
ceiling(max(catchNL$lon)+0.5)),
ylim=c(floor(min(catchNL$lat)-0.5),
ceiling(max(catchNL$lat)+0.5)))+
geom_scatterpie(aes(x=lon, y=lat, group = HaulNumber, r = 0.1),
data = catchNL, cols = speciesNames)
save(CatchAll, HaulAll, BiologyAll, CruiseAll, DataAll, file = file.path(rawDataPath,paste0(surveyYear,"_HERAS",".Rdata")))
}
......
......@@ -37,10 +37,10 @@ resultsPath <- file.path(".","results")
source(file.path(functionPath,"fill_missing_species_v3.R"))
source(file.path(functionPath,"fill_missing_mat_v1.R"))
source(file.path(functionPath,"build_stox_project_tree.R"))
source(file.path(functionPath,"compute_NSAS_WBSS_strata.R"))
source(file.path(functionPath,"compute_NSAS_WBSS_strata_mat.R"))
load(file.path(dataPath,'split_prop_NO.Rdata'))
project_mapping <- read.csv(file.path(dataPath,"project_mapping.csv"))
project_mapping <- read.csv(file.path(dataPath,"2_project_mapping_baseline.csv"))
surveyYearMat <- unique(project_mapping$year)
......@@ -48,7 +48,7 @@ nboot <- 500
# parameters HERAS index
minYear <- 1989
maxYear <- 2019
maxYear <- 2020
yearVec <- minYear:maxYear
nYears <- length(yearVec)
minAge <- 0
......@@ -56,6 +56,8 @@ maxAge <- 9
ageVec <- minAge:maxAge
nAges <- length(ageVec)
nIter <- 1
strataAll <- c(31,41,21,71,151,42,152,111,121,91,51,131,61,81,101,11,141)
nStrata <- length(strataAll)
######################################################################
######################################################################
......@@ -66,11 +68,11 @@ nIter <- 1
######################################################################
# stock object
FLR_temp <- FLQuant( array(-1,dim = c(nAges,nYears,1,1,1,nIter)),
FLR_temp <- FLQuant( array(0,dim = c(nAges,nYears,1,1,1,nIter)),
dimnames=list( age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season="all",
season='all',
area='all'),
units='1e6')
......@@ -91,23 +93,23 @@ WBSS <- FLStock(name = 'WBSS',
######################################################################
# EU NO index object
FLR_HERAS_EU_NO <- FLQuant( array(-1,dim = c(nAges,nYears,1,1,2,nIter)),
dimnames=list(age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season="all",
area=c('EU','NO')),
units='1e6')
FLR_HERAS_EU_NO <- FLIndex(name = 'HERAS_EU_NO',
desc = 'HERAS index from EU and NO StoX projects',
index=FLR_HERAS_EU_NO)
FLR_HERAS_strata <- FLQuant( array(0,dim = c(nAges,nYears,1,3,nStrata,nIter)),
dimnames=list(age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season=c('all','IMM','MAT'),
area=ac(strataAll)),
units='1e6')
FLR_HERAS_EU_NO <- FLIndex(name = 'HERAS_strata',
desc = 'HERAS index disagregated by strata',
index=FLR_HERAS_strata)
# total index object
FLR_HERAS <- FLQuant( array(-1,dim = c(nAges,nYears,1,1,1,nIter)),
FLR_HERAS <- FLQuant( array(0,dim = c(nAges,nYears,1,3,1,nIter)),
dimnames=list(age=as.character(ageVec),
year=as.character(yearVec),
unit="1e6",
season="all",
unit=c("1e6"),
season=c('all','IMM','MAT'),
area='all'),
units='1e6')
FLR_HERAS <- FLIndex( name = 'HERAS_all',
......@@ -186,30 +188,34 @@ for(idxProject in 1:dim(project_mapping)[1]){
endTabNO <- currentBaseLineReport$outputData$FillMissingData
endTabNO[is.na(endTabNO)] <- '-' # replace NA with '-'
#################################################################
## create index and assign values
#################################################################
outIndex <- compute_NSAS_WBSS_strata( endTabEU,
endTabNO,
split_prop,
ageVec,
surveyYear)
# assign numbers at age to FLR objects
HERAS.NSAS$HERAS_EU_NO@index[,ac(surveyYear),,,'NO'] <- colSums(outIndex$NSAS_NO)*1e-6
HERAS.WBSS$HERAS_EU_NO@index[,ac(surveyYear),,,'NO'] <- colSums(outIndex$WBSS_NO)*1e-6
HERAS.NSAS$HERAS_EU_NO@index[,ac(surveyYear),,,'EU'] <- colSums(outIndex$NSAS_EU)*1e-6
HERAS.WBSS$HERAS_EU_NO@index[,ac(surveyYear),,,'EU'] <- colSums(outIndex$WBSS_EU)*1e-6
outIndex <- compute_NSAS_WBSS_strata_mat( endTabEU,
endTabNO,
split_prop,
ageVec,
surveyYear,
strataAll)
for(idxStrata in strataAll){
for(idxMAT in c('all','IMM','MAT')){
HERAS.NSAS$HERAS_strata@index[,ac(surveyYear),,idxMAT,ac(idxStrata)] <-
outIndex$NSAS[ac(idxStrata),,idxMAT]*1e-6
HERAS.WBSS$HERAS_strata@index[,ac(surveyYear),,idxMAT,ac(idxStrata)] <-
outIndex$WBSS[ac(idxStrata),,idxMAT]*1e-6
}
}
# fill in stock weight at age and maturity
NSAS@stock.wt[,ac(surveyYear)] <- outIndex$wa$NSAS
NSAS@mat[,ac(surveyYear)] <- outIndex$propM$NSAS
NSAS@stock.n[,ac(surveyYear)] <- areaSums(HERAS.NSAS$HERAS_EU_NO@index[,ac(surveyYear)])
NSAS@stock.n[,ac(surveyYear)] <- areaSums(HERAS.NSAS$HERAS_strata@index[,ac(surveyYear),,'all'])
WBSS@stock.wt[,ac(surveyYear)] <- outIndex$wa$WBSS
WBSS@mat[,ac(surveyYear)] <- outIndex$propM$WBSS
WBSS@stock.n[,ac(surveyYear)] <- areaSums(HERAS.WBSS$HERAS_EU_NO@index[,ac(surveyYear)])
WBSS@stock.n[,ac(surveyYear)] <- areaSums(HERAS.WBSS$HERAS_strata@index[,ac(surveyYear),,'all'])
# fill in SSB
NSAS@stock[,ac(surveyYear)] <- sum(NSAS@stock.n[,ac(surveyYear)]*NSAS@stock.wt[,ac(surveyYear)]*NSAS@mat[,ac(surveyYear)])
......@@ -223,8 +229,8 @@ setwd('..')
# combine NO and EU components
HERAS.NSAS$HERAS_all@index[] <- areaSums(HERAS.NSAS$HERAS_EU_NO@index)
HERAS.WBSS$HERAS_all@index[] <- areaSums(HERAS.WBSS$HERAS_EU_NO@index)
HERAS.NSAS$HERAS_all@index[] <- areaSums(unitSums(HERAS.NSAS$HERAS_strata@index))
HERAS.WBSS$HERAS_all@index[] <- areaSums(unitSums(HERAS.WBSS$HERAS_strata@index))
HERAS_baseline <- list(HERAS.NSAS=HERAS.NSAS,
HERAS.WBSS=HERAS.WBSS,
......
# script that derives the HERAS index directly from the StoX projects.
rm(list=ls())
library(Rstox)
library(FLFleet)
library(FLCore)
library(ggplotFL)
library(mgcv)
library(tidyverse)
library(stats)
library(XML)
library(tidyr)
#Set up directories
#path <- 'J:/git/HERAS/'
path <- 'J:/git/heras_index_kbwot'
try(setwd(path),silent=TRUE)
mainPath <- file.path(".")
dataPath <- file.path(".","data")
stoxPath <- file.path(".","StoX")
outPath <- file.path(".","output")
functionPath <- file.path(".","functions")
resultsPath <- file.path(".","results")
######################################################################
######################################################################
######################################################################
#### Define parameters and load functions and data
######################################################################
######################################################################
######################################################################
source(file.path(functionPath,"fill_missing_species_v3.R"))
source(file.path(functionPath,"fill_missing_mat_v1.R"))
source(file.path(functionPath,"build_stox_project_tree.R"))
source(file.path(functionPath,"compute_NSAS_WBSS_strata_mat.R"))
load(file.path(dataPath,'split_prop_NO.Rdata'))
project_mapping <- read.csv(file.path(dataPath,"project_mapping.csv"))
surveyYearMat <- unique(project_mapping$year)
nboot <- 500
# parameters HERAS index
minYear <- 1989
maxYear <- 2020
yearVec <- minYear:maxYear
nYears <- length(yearVec)
minAge <- 0
maxAge <- 9
ageVec <- minAge:maxAge
nAges <- length(ageVec)
nIter <- 1
strataAll <- c(31,41,21,71,151,42,152,111,121,91,51,131,61,81,101,11,141)
nStrata <- length(strataAll)
######################################################################
######################################################################
######################################################################
#### # initialise FLStock objects
######################################################################
######################################################################
######################################################################
# stock object
FLR_temp <- FLQuant( array(0,dim = c(nAges,nYears,1,1,1,nIter)),
dimnames=list( age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season='all',
area='all'),
units='1e6')
NSAS <- FLStock(name = 'NSAS',
desc = 'NSAS',
FLR_temp)
WBSS <- FLStock(name = 'WBSS',
desc = 'WBSS',
FLR_temp)
######################################################################
######################################################################
######################################################################
#### # initialise FLIndex objects
######################################################################
######################################################################
######################################################################
# EU NO index object
FLR_HERAS_strata <- FLQuant( array(0,dim = c(nAges,nYears,1,3,nStrata,nIter)),
dimnames=list(age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season=c('all','IMM','MAT'),
area=ac(strataAll)),
units='1e6')
FLR_HERAS_EU_NO <- FLIndex(name = 'HERAS_strata',
desc = 'HERAS index disagregated by strata',
index=FLR_HERAS_strata)
# total index object
FLR_HERAS <- FLQuant( array(0,dim = c(nAges,nYears,1,3,1,nIter)),
dimnames=list(age=as.character(ageVec),
year=as.character(yearVec),
unit=c("1e6"),
season=c('all','IMM','MAT'),
area='all'),
units='1e6')
FLR_HERAS <- FLIndex( name = 'HERAS_all',
desc = 'HERAS index (total)',
index=FLR_HERAS)
HERAS.NSAS <- FLIndices(FLR_HERAS_EU_NO, FLR_HERAS)
HERAS.WBSS <- FLIndices(FLR_HERAS_EU_NO, FLR_HERAS)
######################################################################
######################################################################
######################################################################
#### Compute NSAS/WBSS indices
######################################################################
######################################################################
######################################################################
dir.create("StoX",showWarnings = FALSE)
setwd('StoX')
for(idxProject in 1:dim(project_mapping)[1]){
surveyYear <- project_mapping$year[idxProject]
print(project_mapping$project_EU[idxProject])
dir.create(as.character(surveyYear),showWarnings = FALSE)
setwd(as.character(surveyYear))
#################################################################
## create project
#################################################################
# EU project
build_stox_project_tree(project_mapping$project_EU[idxProject],
surveyYear)
runBaseline(projectName = file.path(getwd(),
project_mapping$project_EU[idxProject]),
save = TRUE,
exportCSV = TRUE,
modelType = c('baseline'))
runBaseline(projectName = file.path(getwd(),
project_mapping$project_EU[idxProject]),
save = TRUE,
exportCSV = TRUE,
modelType = c('baseline-report'))
currentBaseLineReport <- getBaseline(file.path(getwd(),
project_mapping$project_EU[idxProject]),
save = FALSE,
exportCSV = FALSE,
modelType = c('baseline-report'),reset = FALSE)
endTabEU <- currentBaseLineReport$outputData$FillMissingData
endTabEU[is.na(endTabEU)] <- '-' # replace NA with '-'
# NO project
build_stox_project_tree(project_mapping$project_NO[idxProject],
surveyYear)
runBaseline(projectName = file.path(getwd(),
project_mapping$project_NO[idxProject]),
save = TRUE,
exportCSV = TRUE,
modelType = c('baseline'))
runBaseline(projectName = file.path(getwd(),
project_mapping$project_NO[idxProject]),
save = TRUE,
exportCSV = TRUE,
modelType = c('baseline-report'))
currentBaseLineReport <- getBaseline(file.path(getwd(),
project_mapping$project_NO[idxProject]),
save = FALSE,
exportCSV = FALSE,
modelType = c('baseline-report'),reset = FALSE)