From 5370460c4b809e3daa577d48c71cafb87df0fac2 Mon Sep 17 00:00:00 2001
From: nielshintzen <nhintzen@pelagicfish.eu>
Date: Thu, 25 Apr 2024 10:40:28 +0200
Subject: [PATCH] adding functionality for set.pars to simulate function

---
 R/SAM2FLR.R | 18 +++++++++++++++++-
 R/methods.R | 14 +++++++++-----
 2 files changed, 26 insertions(+), 6 deletions(-)

diff --git a/R/SAM2FLR.R b/R/SAM2FLR.R
index fbb02dd..45b30bc 100755
--- a/R/SAM2FLR.R
+++ b/R/SAM2FLR.R
@@ -151,7 +151,7 @@ SAM2FLR <- function(fit,ctrl){
   res@range   <- ctrl@range
 return(res)}
 
-SIM2FLR <- function(sim,fit,ctrl){
+SIM2FLR <- function(sim,fit,ctrl,set.pars=NULL){
   resList       <-list()
   for(i in 1:nrow(sim)){
     #- Create new FLSAM object and fill with ctrl elements
@@ -167,6 +167,22 @@ SIM2FLR <- function(sim,fit,ctrl){
     #- get parameter names
     pars2match <- unique(dimnames(sim)[[2]])
     keepPars    <- names(fit$pl)[which(names(fit$pl) %in% pars2match)]
+
+    #- Fix for set.pars
+    if(!is.null(set.pars)){
+      matchNames <- matrix(c("logN.vars","logSdLogN",
+                             "logP.vars","logSdLogP",
+                             "catchabilities","logFpar",
+                             "power.law.exps","logQpow",
+                             "f.vars","logSdLogFsta",
+                             "obs.vars","logSdLogObs"),ncol=2,byrow=T,dimnames=list(1:6,c("FLSAM","SAM")))
+      if(any(!set.pars$par %in% matchNames[,"FLSAM"]))
+          warning(paste(set.pars$par[which(!set.pars$par %in% matchNames[,"FLSAM"])],"cannot be set"))
+      set.parsSAM <- merge(set.pars,matchNames,by.x="par",by.y="FLSAM")
+      for(j in 1:nrow(set.pars))
+        fit$pl[[set.parsSAM$SAM[j]]] <- fit$pl[[set.parsSAM$SAM[j]]][-(set.parsSAM$number[j]+1)]
+    }
+
     pars       <- ac(rownames(as.data.frame(unlist(fit$pl))))
     idx       <- character()
     for(iK in keepPars){
diff --git a/R/methods.R b/R/methods.R
index 295c374..2e457a0 100644
--- a/R/methods.R
+++ b/R/methods.R
@@ -382,12 +382,15 @@ setMethod("coefficients",signature(object="FLSAMs"),
         }
 )
 
-setGeneric("simulate", function(x,y,z,n) standardGeneric("simulate"))
+setGeneric("simulate", function(x,y,z,n,set.pars) standardGeneric("simulate"))
 
 setMethod("simulate",signature(x="FLStock",y="FLIndices",z="FLSAM.control",
-  n='numeric'),
-          function(x,y,z,n=100){
-            fit   <- FLSAM(x,y,z,return.fit=T)
+  n='numeric',set.pars="data.frame"),
+          function(x,y,z,n=100,set.pars=NULL){
+            if(!is.null(set.pars))
+              fit   <- FLSAM(x,y,z,return.fit=T,set.pars=set.pars)
+            if(is.null(set.pars))
+              fit   <- FLSAM(x,y,z,return.fit=T)
             sdrep <- sdreport(fit$obj,getJointPrecision=T)
             sigma <- as.matrix(solve(sdrep$jointPrecision))
             mu    <- c(sdrep$par.fixed,sdrep$par.random)
@@ -395,7 +398,8 @@ setMethod("simulate",signature(x="FLStock",y="FLIndices",z="FLSAM.control",
             colnames(sim) <- rownames(sigma)
             rownames(sim) <- 1:n
             
-            sim   <- SIM2FLR(sim,fit,z)
+            sim   <- SIM2FLR(sim,fit,z,set.pars)
+
           return(sim)})