jags.model <- function(file, data=sys.frame(sys.parent()), inits,
                       n.chains = 1, n.adapt=1000, nchain)
{

    if (missing(file)) {
        stop("Model file name missing")
    }
    if (!file.exists(file)) {
        stop(paste("Model file \"", file, "\" not found", sep=""))
    }
    if (!missing(nchain)) {
        warning("Argument nchain in jags.model is deprecated. Use n.chains.")
        if (missing(n.chains)) {
            n.chains = nchain
        }
    }
    
    p <- .Call("make_console", PACKAGE="CoGAPS") 
    .Call("check_model", p, file, PACKAGE="CoGAPS")

    varnames <- .Call("get_variable_names", p, PACKAGE="CoGAPS")
    if (is.environment(data)) {
        ##Get a list of numeric objects from the supplied environment
        data <- mget(varnames, envir=data, mode="numeric",
                     ifnotfound=list(NULL))
        ##Strip null entries
        data <- data[!sapply(data, is.null)]
    }
    else if (is.list(data)) {
        v <- names(data)
        if (is.null(v)) {
            stop("data must be a named list")
        }
        if (any(nchar(v)==0)) {
            stop("unnamed variables in data list")
        }
        if (any(duplicated(v))) {
            stop("Duplicated names in data list: ",
                 paste(v[duplicated(v)], collapse=" "))
        }
        relevant.variables <- v %in% varnames
        data <- data[relevant.variables]
        unused.variables <- setdiff(v, varnames)
        for (i in seq(along=unused.variables)) {
            warning("Unused variable \"", unused.variables[i], "\" in data")
        }
    }
    else {
        stop("data must be a list or environment")
    }
    
    .Call("compile", p, data, as.integer(n.chains), TRUE, PACKAGE="CoGAPS")

### Setting initial values

    if (!missing(inits)) {

        checkParameters <- function(inits) {
            if(!is.list(inits))
                return (FALSE)

            inames <- names(inits)
            if (is.null(inames) || any(nchar(inames) == 0))
                return (FALSE)

            if (any(duplicated(inames)))
                return (FALSE)
            
            if (any(inames==".RNG.name")) {
                rngname <- inits[[".RNG.name"]]
                if (!is.character(rngname) || length(rngname) != 1)
                    return (FALSE)
                inits[[".RNG.name"]] <- NULL
            }

            if (!all(sapply(inits, is.numeric)))
                return (FALSE)
            
            return (TRUE)
        }
        
        setParameters <- function(inits, chain) {
            if (!is.null(inits[[".RNG.name"]])) {
                .Call("set_rng_name", p, inits[[".RNG.name"]],
                      as.integer(chain), PACKAGE="CoGAPS")
                inits[[".RNG.name"]] <- NULL
            }
            .Call("set_parameters", p, inits, as.integer(chain),
                  PACKAGE="CoGAPS")
        }
        
        init.values <- vector("list", n.chains)
        
        if (is.function(inits)) {
            if (any(names(formals(inits)) == "chain")) {
                for (i in 1:n.chains) {
                    init.values[[i]] <- inits(chain=i)
                }
            }
            else {
                for (i in 1:n.chains) {
                    init.values[[i]] <- inits()
                }
            }
        }
        else if (is.list(inits)) {

            if (checkParameters(inits)) {
                ## Replicate initial values for all chains
                for (i in 1:n.chains) {
                    init.values[[i]] <- inits
                }
            }
            else {
                if (length(inits) != n.chains) {
                    stop("Length mismatch in inits")
                }
                init.values <- inits
            }
        }
            
        for (i in 1:n.chains) {
            if (!checkParameters(init.values[[i]])) {
                stop("Invalid parameters for chain ", i)
            }
            setParameters(init.values[[i]], i)
        }
    }

    .Call("initialize", p, PACKAGE="CoGAPS")

    model.state <- .Call("get_state", p, PACKAGE="CoGAPS")
    model.data <- .Call("get_data", p, PACKAGE="CoGAPS")
    model.code <- readLines(file, warn=FALSE)
    model <- list("ptr" = function() {p},
                  "data" = function() {model.data},
                  "model" = function() {model.code},
                  "state" = function(internal=FALSE)
                  {
                      if(!internal) {
                          for(i in 1:n.chains) {
                              model.state[[i]][[".RNG.state"]] <- NULL
                              model.state[[i]][[".RNG.name"]] <- NULL
                          }
                      }
                      return(model.state)
                  },
                  "nchain" = function()
                  {
                      .Call("get_nchain", p, PACKAGE="CoGAPS")
                  },
                  "iter" = function()
                  {
                      .Call("get_iter", p, PACKAGE="CoGAPS")
                  },
                  "sync" = function() {
                      
                      model.state <<- .Call("get_state", p, PACKAGE="CoGAPS")
                  },
                  "recompile" = function() {
                      ## Clear the console
                      .Call("clear_console", p, PACKAGE="CoGAPS")
                      p <<- .Call("make_console", PACKAGE="CoGAPS")
                      ## Write the model to a temporary file so we can re-read it
                      mf <- tempfile()
                      writeLines(model.code, mf)
                      .Call("check_model", p, mf, PACKAGE="CoGAPS")
                      unlink(mf)
                      ## Re-compile
                      .Call("compile", p, data, n.chains, FALSE, PACKAGE="CoGAPS")
                      ## Re-initialize
                      if (!is.null(model.state)) {
                          if (length(model.state) != n.chains) {
                              stop("Incorrect number of chains in saved state")
                          }
                          for (i in 1:n.chains) {
                              statei <- model.state[[i]]
                              rng <- statei[[".RNG.name"]]
                              if (!is.null(rng)) {
                                  .Call("set_rng_name", p, rng, i, PACKAGE="CoGAPS")
                                  statei[[".RNG.name"]] <- NULL
                              }
                              .Call("set_parameters", p, statei, i, PACKAGE="CoGAPS")
                          }
                          .Call("initialize", p, PACKAGE="CoGAPS")
                          ## Redo adaptation
                          adapting <- .Call("is_adapting", p, PACKAGE="CoGAPS")
                          if(n.adapt > 0 && adapting) {
                              cat("Adapting\n")
                              .Call("update", p, n.adapt, PACKAGE="CoGAPS")
                              if (!.Call("adapt_off", p, PACKAGE="CoGAPS")) {
                                  warning("Adaptation incomplete");
                              }
                          }
                          model.state <<- .Call("get_state", p, PACKAGE="CoGAPS")
                      }
                      invisible(NULL)
                  })
    class(model) <- "jags"

    if (n.adapt > 0) {
        adapt(model, n.adapt)
    }
    return(model)
}

parse.varname <- function(varname) {

  ## Try to parse string of form "a" or "a[n,p:q,r]" where "a" is a
  ## variable name and n,p,q,r are integers

  v <- try(parse(text=varname, n=1), silent=TRUE)
  if (!is.expression(v) || length(v) != 1)
    return(NULL)

  v <- v[[1]]
  if (is.name(v)) {
    ##Full node array requested
    return(list(name=deparse(v)))
  }
  else if (is.call(v) && identical(deparse(v[[1]]), "[") && length(v) > 2) {
    ##Subset requested
    ndim <- length(v) - 2
    lower <- upper <- numeric(ndim)
    if (any(nchar(sapply(v, deparse)) == 0)) {
      ##We have to catch empty indices here or they will cause trouble
      ##below
      return(NULL)
    }
    for (i in 1:ndim) {
      index <- v[[i+2]]
      if (is.numeric(index)) {
        ##Single index
        lower[i] <- upper[i] <- index
      }
      else if (is.call(index) && length(index) == 3 &&
               identical(deparse(index[[1]]), ":") &&
               is.numeric(index[[2]]) && is.numeric(index[[3]]))
        {
          ##Index range
          lower[i] <- index[[2]]
          upper[i] <- index[[3]]
        }
      else return(NULL)
    }
    if (any(upper < lower))
      return (NULL)
    return(list(name = deparse(v[[2]]), lower=lower, upper=upper))
  }
  return(NULL)
}

parse.varnames <- function(varnames)
{
  names <- character(length(varnames))
  lower <- upper <- vector("list", length(varnames))
  for (i in seq(along=varnames)) {
    y <- parse.varname(varnames[i])
    if (is.null(y)) {
      stop(paste("Invalid variable subset", varnames[i]))
    }
    names[i] <- y$name
    if (!is.null(y$lower)) {
      lower[[i]] <- y$lower
    }
    if (!is.null(y$upper)) {
      upper[[i]] <- y$upper
    }
  }
  return(list(names=names, lower=lower, upper=upper))
}


jags.samples <-
  function(model, variable.names, n.iter, thin=1, type="trace", ...)
{
    if (class(model) != "jags")
      stop("Invalid JAGS model")

    if (!is.character(variable.names) || length(variable.names) == 0)
      stop("variable.names must be a character vector")
     
    if (!is.numeric(n.iter) || length(n.iter) != 1 || n.iter <= 0)
      stop("n.iter must be a positive integer")
    if (!is.character(type))
      stop("type must be a character vector")

    pn <- parse.varnames(variable.names)
    .Call("set_monitors", model$ptr(), pn$names, pn$lower, pn$upper,
          as.integer(thin), type, PACKAGE="CoGAPS")
    update(model, n.iter, ...)
    ans <- .Call("get_monitored_values", model$ptr(), type, PACKAGE="CoGAPS")
    for (i in seq(along=variable.names)) {
      .Call("clear_monitor", model$ptr(), pn$names[i], pn$lower[[i]],
            pn$upper[[i]], type, PACKAGE="CoGAPS")
    }
    return(ans)
}

nchain <- function(model)
{
    if (!inherits(model, "jags"))
      stop("Invalid JAGS model object in nchain")
    
    .Call("get_nchain", model$ptr(), PACKAGE="CoGAPS")
}

load.module <- function(name, path, quiet=FALSE)
{
    if (name %in% list.modules()) {
        ## This is a stop-gap measure as JAGS 2.1.0 does allow you
        ## to load the same module twice. This should be fixed in
        ## later versions.
        return(invisible()) #Module already loaded
    }
    
    if (missing(path)) {
        path = getOption("jags.moddir")
        if (is.null(path)) {
            stop("option jags.moddir is not set")
        }
    }
    if (!is.character(path) || length(path) != 1)
        stop("invalid path")
    if (!is.character(name) || length(name) != 1)
        stop("invalid name")

    file <- file.path(path, paste(name, .Platform$dynlib.ext, sep=""))
    if (!file.exists(file)) {
        stop("File not found: ", file)
    }
    if (!isDLLLoaded(file)) {
        ## We must avoid calling dyn.load twice on the same DLL This
        ## may result in the DLL being unloaded and then reloaded,
        ## which will invalidate pointers to the distributions,
        ## functions and factories in the module.
        dyn.load(file)
    }
    ok <- .Call("load_module", name, PACKAGE="CoGAPS")
    if (!ok) {
        stop("module", name, "not found\n", sep=" ")
    }
    else if (!quiet) {
        cat("module", name, "loaded\n", sep=" ")
    }
    invisible()
}

list.modules <- function()
{
    .Call("get_modules", PACKAGE="CoGAPS");
}

isDLLLoaded <- function(file)
{
    dll.list <- getLoadedDLLs()
    for (i in seq(along=dll.list)) {
        if (dll.list[[i]]["path"][1] == file)
            return(TRUE)
    }
    return(FALSE)
}