library(stockassessment)
library(parallel)

thisrunwithout <- function (fit, year = NULL, fleet = NULL, map = fit$obj$env$map, ...) {
    data <- reduce(fit$data, year = year, fleet = fleet, conf = fit$conf)
    conf <- attr(data, "conf")
    fakefile <- file()
    sink(fakefile)
    saveConf(conf, file = "")
    sink()
    conf <- loadConf(data, fakefile, patch = TRUE)
    close(fakefile)
    par <- defpar(data, conf)
    ret <- sam.fit(data, conf, par, rm.unidentified = TRUE, map = map, lower = fit$low, upper = fit$hig, ...)
    return(ret)
}


thisretro <- function (fit, year = NULL, ncores = detectCores(), ...) 
{
    data <- fit$data
    y <- fit$data$aux[, "year"]
    f <- fit$data$aux[, "fleet"]
    suf <- sort(unique(f))
    maxy <- sapply(suf, function(ff) max(y[f == ff]))
    if (length(year) == 1) {
        mat <- sapply(suf, function(ff) {
            my <- maxy[ff]
            my:(my - year + 1)
        })
        if (year == 1) 
            mat <- matrix(mat, nrow = 1)
    }
    if (is.vector(year) & length(year) > 1) {
        mat <- sapply(suf, function(ff) year)
    }
    if (is.matrix(year)) {
        mat <- year
    }
    if (nrow(mat) > length(unique(y))) 
        stop("The number of retro runs exceeds number of years")
    if (ncol(mat) != length(suf)) 
        stop("Number of retro fleets does not match")
    setup <- lapply(1:nrow(mat), function(i) do.call(rbind, lapply(suf, 
        function(ff) if (mat[i, ff] <= maxy[ff]) {
            cbind(mat[i, ff]:maxy[ff], ff)
        })))
    if (ncores > 1) {
        cl <- makeCluster(ncores)
        on.exit(stopCluster(cl))
        clusterExport(cl, varlist = "fit", envir = environment())
        clusterExport(cl, varlist = "thisrunwithout", envir = environment())
        lib.ver <- dirname(path.package("stockassessment"))
        clusterExport(cl, varlist = "lib.ver", envir = environment())
        clusterEvalQ(cl, {
            library(stockassessment, lib.loc = lib.ver)
        })
        runs <- parLapply(cl, setup, function(s) thisrunwithout(fit, 
            year = s[, 1], fleet = s[, 2], ...))
    }
    else {
        runs <- lapply(setup, function(s) thisrunwithout(fit, year = s[, 
            1], fleet = s[, 2], ...))
    }
    converg <- unlist(lapply(runs, function(x) x$opt$conv))
    if (any(converg != 0)) 
        warning(paste0("retro run(s) ", paste0(which(converg != 
            0), collapse = ","), " did not converge."))
    attr(runs, "fit") <- fit
    class(runs) <- "samset"
    runs
}







load("run/model.RData")
RETRO<-thisretro(fit, year=5)
save(RETRO, file="run/retro.RData")
