### A collection of scripts for dealing with ECD model structures.


##Creates a 'daCapo' object which consists of a collection of
##mcmc.list objects organized by model and variable.
"daCapo" <-
function(chains, window, studentModelList = listStudentModels(),
         linkModelList = listLinkModels(),trimRegExp=NULL,
         maxparams=99) {
  ## Argument defaults
  if (missing(chains)) {
    chains <- getChains()[getCycles()[,"chain_last_cycle"]>0]
  }
  if (missing(window)) {
    maxcycle <- min(getCycles()[chains,"chain_last_cycle"])
    burnin <- 1001
    if (maxcycle < 3000) burnin <- 501
    if (maxcycle < 600) burnin <- 1
    window <- c(burnin,maxcycle)
  }
  ## If trimRegExp is non-null, then remove matching expressions from
  ## student model list and link model list
  if (!is.null(trimRegExp)) {
    trimindex <- grep(trimRegExp,studentModelList)
    if (length(trimindex)>0) {
      studentModelList <- studentModelList [-trimindex]
    }
    trimindex <- grep(trimRegExp,linkModelList)
    if (length(trimindex)>0) {
      linkModelList <- linkModelList [-trimindex]
    }
  }
  studentModels.pnames <- sapply(studentModelList,listParameterNames)
  linkModels.pnames <- sapply(linkModelList,listParameterNames)
  ###Load parameters
  studentModels.parameters <-
    lapply(studentModelList, getAllParameters,
           iterations=window,chain=chains,maxparams=maxparams)
  names(studentModels.parameters) <- studentModelList
  linkModels.parameters <-
    lapply(linkModelList, getAllParameters,
           iterations=window,chain=chains,maxparams=maxparams)
  names(linkModels.parameters) <- linkModelList
  result <- list(studentModels=studentModels.parameters,
                 linkModels=linkModels.parameters,
                 chains,window)
  class(result) <- "daCapo"
  result
}

  

## Applies a function recursively to the parts of the daCapo object
## and prints the result. 
## The arguments whichStudentModels and whichLinkModels can be used
## to print out only a subset of the models.
"daCapoApplyPrint" <-
function(daCapo.obj, FUN,...,
         whichStudentModels=1:length(daCapo.obj$studentModels),
         whichLinkModels=1:length(daCapo.obj$linkModels)) {
  FUN <- match.fun(FUN)
  for (model in names(daCapo.obj$studentModels[whichStudentModels])) {
    cat ("\nStudent Model ",model,"\n")
    paramlist <- daCapo.obj$studentModels[[model]]
    for (pname in names(paramlist)) {
      cat("Values for parameter ",pname,"\n")
      param <- paramlist[[pname]]
      print(try(FUN(param,...)))
      cat("\n")
    }
    cat("\n")
  }

  for (model in names(daCapo.obj$linkModels[whichLinkModels])) {
    cat ("\nLink Model ",model,"\n")
    paramlist <- daCapo.obj$linkModels[[model]]
    for (pname in names(paramlist)) {
      cat("Values for parameter ",pname,"\n")
      param <- paramlist[[pname]]
      print(try(FUN(param,...)))
      cat("\n")
    }
    cat("\n")
  }

}



## Applies a function recursively to the parts of the daCapo object
## and assumes each result returns a plot. It adds a sub argument to
## the plot based on the model.
## The arguments whichStudentModels and whichLinkModels can be used
## to print out only a subset of the models.
## If prefix is supplied, each model is plotted to a separate PDF file,
## If ask is true, plots are made with
"daCapoApplyPlot" <-
function(daCapo.obj, FUN,...,
         whichStudentModels=1:length(daCapo.obj$studentModels),
         whichLinkModels=1:length(daCapo.obj$linkModels),
         prefix="",ask=missing(prefix)) {
  FUN <- match.fun(FUN)
  onscreen <- missing(prefix)
  first <- TRUE
  for (model in names(daCapo.obj$studentModels[whichStudentModels])) {
    if (!onscreen) {
      pdf(paste(prefix,model,".pdf",sep=""))
    }
    if (!first && ask) {
      readline("Hit return for next model: ")
    }
    first <- TRUE
    paramlist <- daCapo.obj$studentModels[[model]]
    for (pname in names(paramlist)) {
      param <- paramlist[[pname]]
      if (!first && ask) {
        readline("Hit return for next variable: ")
      }
      first <- FALSE
      try(FUN(param,...,ask=ask,
              sub=paste("Model: ",model,"Par: ",getObservableVar(pname),sep=" ")))
    }
    if (!onscreen) {
      dev.off()
    }
  }

  for (model in names(daCapo.obj$linkModels[whichLinkModels])) {
    if (!onscreen) {
      pdf(paste(prefix,model,".pdf",sep=""))
    }
    if (!first && ask) {
      readline("Hit return for next model: ")
    }
    first <- TRUE
    paramlist <- daCapo.obj$linkModels[[model]]
    lmname <- splitlmn(model)
    for (pname in names(paramlist)) {
      param <- paramlist[[pname]]
      if (!first && ask) {
        readline("Hit return for next variable: ")
      }
      first <- FALSE
      try(FUN(param,...,ask=ask,
              sub=paste("Task: ",lmname$task,"Par: ",getObservableVar(pname),sep=" ")))
    }
    if (!onscreen) {
      dev.off()
    }
  }

}