library(stockassessment)
load("run/model.RData")
FC<-list()
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,1,1), label="SQ all years", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,NA,NA,NA), fval=c(NA,0.000001,0.000001, 0.000001), label="SQ then zero", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,NA,NA,NA), fval=c(NA,0.23,0.23, 0.23), label="SQ then Fmsy", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,NA,NA,NA), fval=c(NA,0.60,0.60, 0.60), label="SQ then Fp05", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,NA,NA,NA), fval=c(NA,0.62,0.62, 0.62), label="SQ then Fpa", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,NA,NA,NA), fval=c(NA,0.85,0.85, 0.85), label="SQ then Flim", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)

Ry<-1997:2024     # assessment year-1 
Fsq<-0.685        # Fsq latest value from SAM
yearbase<- 2025   # assessment year-1
 # gives exactly Fsq value
set.seed(12345)
FC[[length(FC)+1]]  <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq, Fsq,Fsq,Fsq),rec.years=Ry, label="Rec 1997-2024, SQ all years", processNoise=FALSE,addTSB=TRUE, year.base=yearbase, savesim=TRUE) 
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.000001, 0.000001), rec.years=Ry, label="Rec 1997-2024, SQ then zero", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.23, 0.23), rec.years=Ry, label="Rec 1997-2024, SQ then Fmsy", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.60, 0.60), rec.years=Ry, label="Rec 1997-2024, SQ then Fpa/Fp05", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.62, 0.62), rec.years=Ry, label="Rec 1997-2024, SQ then Fpa", processNoise=FALSE, addTSB=TRUE, year.base=2023)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.85, 0.85), rec.years=Ry, label="Rec 1997-2024, SQ then Flim", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.128, 0.128), rec.years=Ry, label="Rec 1997-2024, SQ then Fmsy*SSB2025/MSYBtrigger", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)

Ry<-2020:2024     # latest 5 years as recommended from ADGdeep 2022 (wgdeep2023) 
Fsq<-0.685       # Fsq latest value from SAM 
yearbase<- 2025   # assessment year-1
 # gives exactly Fsq value
set.seed(12345)
FC[[length(FC)+1]]  <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq, Fsq,Fsq,Fsq),rec.years=Ry, label="Rec 2020-2024, SQ all years", processNoise=FALSE,addTSB=TRUE, year.base=yearbase, savesim=TRUE) 
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.000001, 0.000001), rec.years=Ry, label="Rec 2020-2024, SQ then zero", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.23, 0.23), rec.years=Ry, label="Rec 2020-2024, SQ then Fmsy", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.60, 0.60), rec.years=Ry, label="Rec 2020-2024, SQ then Fpa/Fp05", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.62, 0.62), rec.years=Ry, label="Rec 2020-2024, SQ then Fpa", processNoise=FALSE, addTSB=TRUE, year.base=2023)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.85, 0.85), rec.years=Ry, label="Rec 2020-2024, SQ then Flim", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA), fval=c(Fsq,Fsq,0.127, 0.127), rec.years=Ry, label="Rec 2020-2024, SQ then Fmsy*SSB2024/MSYBtrigger", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=TRUE)
# Gives Fsq close but not exactly
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,1,1), rec.years=Ry, label="Rec 1996-2020, SQ all years", processNoise=FALSE, addTSB=TRUE, year.base=2021)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,NA,NA), fval=c(NA,NA,0.000001, 0.000001), rec.years=Ry, label="Rec 1996-2020, SQ then zero", processNoise=FALSE, addTSB=TRUE, year.base=2021)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,NA,NA), fval=c(NA,NA,0.23, 0.23), rec.years=Ry, label="Rec 1996-2020, SQ then Fmsy", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,NA,NA), fval=c(NA,NA,0.60, 0.60), rec.years=Ry, label="Rec 1996-2020, SQ then Fp05", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,NA,NA), fval=c(NA,NA,0.62, 0.62), rec.years=Ry, label="Rec 1996-2020, SQ then Fpa", processNoise=FALSE, addTSB=TRUE)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(1,1,NA,NA), fval=c(NA,NA,0.85, 0.85), rec.years=Ry, label="Rec 1996-2020, SQ then Flim", processNoise=FALSE, addTSB=TRUE)


Ry<-2020:2024     # latest 5 years as recommended from ADGdeep 2022 (wgdeep2023) 
Fsq<-0.360        # value for year 2024 
yearbase<- 2025   # assessment year-1
 # gives exactly Fsq value
set.seed(12345)
FC[[length(FC)+1]]  <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq, Fsq,Fsq,Fsq,Fsq,Fsq),rec.years=Ry, label="Rec 2020-2024, SQ all years*", processNoise=FALSE,addTSB=TRUE, year.base=yearbase, savesim=F) 
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.000001, 0.000001,0.000001, 0.000001), rec.years=Ry, label="Rec 2020-2024, SQ then zero*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.23, 0.23,0.23, 0.23), rec.years=Ry, label="Rec 2020-2024, SQ then Fmsy*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.60, 0.60,0.60, 0.60), rec.years=Ry, label="Rec 2020-2024, SQ then Fpa/Fp05*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.62, 0.62), rec.years=Ry, label="Rec 2020-2024, SQ then Fpa", processNoise=FALSE, addTSB=TRUE, year.base=yearbase)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.85, 0.85,0.85, 0.85), rec.years=Ry, label="Rec 2020-2024, SQ then Flim*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.127, 0.127,0.127, 0.127), rec.years=Ry, label="Rec 2020-2024, SQ then Fmsy*SSB2025/MSYBtrigger*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)


## WGDEEP 2026 test
Ry<-2015:2024     # latest 10 years  
Fsq<-0.37        # latest 5 years 
yearbase<- 2025   # assessment year-1
 # gives exactly Fsq value
set.seed(12345)
FC[[length(FC)+1]]  <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq, Fsq,Fsq,Fsq,Fsq,Fsq),rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ all years*", processNoise=FALSE,addTSB=TRUE, year.base=yearbase, savesim=F) 
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.000001, 0.000001,0.000001, 0.000001), rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ then zero*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.23, 0.23,0.23, 0.23), rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ then Fmsy*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.60, 0.60,0.60, 0.60), rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ then Fpa/Fp05*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
#set.seed(12345)
#FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.62, 0.62), rec.years=Ry, label="Rec 2015-2024, SQ then Fpa", processNoise=FALSE, addTSB=TRUE, year.base=yearbase)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.85, 0.85,0.85, 0.85), rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ then Flim*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)
set.seed(12345)
FC[[length(FC)+1]] <- forecast(fit, fscale=c(NA,NA,NA,NA,NA,NA), fval=c(Fsq,Fsq,0.127, 0.127,0.127, 0.127), rec.years=Ry, label="Rec 2015-2024, F=0.37, SQ then Fmsy*SSB2025/MSYBtrigger*", processNoise=FALSE, addTSB=TRUE, year.base=yearbase, savesim=F)






save(FC, file="run/forecast.RData")

