/*********************************************************
** Copyright (c) 2005
** University of Washington
** Licensed under the terms set forth by University of
** Washington. If you did not sign such a license, you
** are using this software/code illegally and you do not
** have permission to use, modify, or redistribute
** this or any files in this software package.
**
** File: RInterfaces.cpp  $Revision: 274 $
**
**********************************************************/
#include <iostream>
#include "RInterfaces.h"
#include "OdpThread.h"
#include "TimecourseThread.h"
#include "TimecourseNoIntercept.h"
#include "TimecourseIntercept.h"
#include "ProgressThreadDef.h"
#include "ProgressThread.h"
#include "RInterfaces.h"

void CheckLength(long actLength, const char* argName, long expLength)
{
  if(actLength != expLength)
    {
      Rf_error("length of '%s' is %d; expected %d", argName, actLength, expLength);
    }
}

#ifdef __cplusplus
extern "C" {
#endif
  /* Internal function headers */
SEXP GetODPResultList (COdpThread *pThread);
SEXP GetTimecourseResultList(CProgressThread *pThread);

  /*
   * Starts an ODP computation in a background thread.
   *
   * Arguments:
   *  xx : Vector of sums of squares
   *  xxk : Matrix of group-specific sums
   *  mu : Matrix of group-specific means
   *  sigma : Vector of alternative standard deviations
   *  sigma0 : Vector of null standard deviations for data
   *  muBS : Variable mu indexed by subset nii
   *  sigmaBS : Variable sigma indexed by subset nii
   *  sigma0BS : Variable sigma0 indexed by subset nii
   *  ID : Model matrix
   *  nGrp : Vector of the number of samples in each group
   *  res : Matrix of residuals
   *  vv : Matrix of indices for bootstrap samples
   *  nii : Vector of subset of genes to use in ODP denominator
   *  m : Number of rows of input data (integer)
   *  n : Number of columns of input data (integer)
   *  k : Number of groups
   *  B : Number of null iterations to compute (integer)
   *  m0 : Number of genes to use in ODP denominator (integer)
   *  startThread : Logical value indicating whether to start computation
   *    as a background thread.
   * Return value:
   *  Integer uniquely identifying the thread that was started.
   */

// gah - change to passing a list instead of separate variables
  SEXP RunODP(SEXP arglist, SEXP startThread)
{
  int nProtected = 0;

  bool bStartThread = LOGICAL_POINTER(startThread)[0];
  int status;

  SEXP tid;
  PROTECT(tid = NEW_INTEGER(1));  nProtected++;

  COdpThread* pOdpThread = NULL;
  try {
    /* Allocate COdpThread object */
    pOdpThread = new COdpThread(arglist);
  } 
  catch(bad_alloc&) {
    pOdpThread = NULL;
  }

  /* Check that object was allocated correctly */
  if(NULL == pOdpThread)
    {
      Rf_error("Could not allocate COdpThread object");
      UNPROTECT(nProtected);
      return R_NilValue;
    }

  if (!bStartThread)
    {
      /* Do computation and return when finished */
      status = pOdpThread->Go(0);
      if(0 != status)
	{
	  /* Computation failed; cleanup and return */
	  delete pOdpThread;
	  Rf_error("ODP computation failed (error number %d)", status);
	  UNPROTECT(nProtected);
	  return R_NilValue;
	}
      else
	{
	  /* Get result data and return */
	  SEXP slist;
	  PROTECT(slist = GetODPResultList(pOdpThread)); nProtected++;
	  delete pOdpThread;
	  UNPROTECT(nProtected);
	  return(slist);
	}
    }
  else
    {
      /* Start background thread and return immediately */
      status = pOdpThread->Start();
      if(0 != status)
	{
	  if(!pOdpThread->GetIsStarted())
	    {
	      /* Thread failed to start */
	      delete pOdpThread;
	      Rf_error("Could not start thread (error number %d)", status);
	      UNPROTECT(nProtected);
	      return R_NilValue;
	    }
	  else
	    {
	      /* Thread started, but something else is wrong */
	      Rf_warning("Thread %d was started, but returned status code %d", 
			 pOdpThread->GetThreadId(), status);
	    }
	}
      
      THREAD_ID_T xtid = pOdpThread->GetThreadId();
      INTEGER_POINTER(tid)[0] = (int)xtid;
      
      UNPROTECT(nProtected);
      return tid;
    }
}

  /* 
   * Gets the result of an ODP background thread that has finished.
   *   Arguments:
   *     tid : ID of thread being queried
   *   Return value: List containing the following:
   *     lr : observed ODP statistics (m-vector)
   *     lr0 : null ODP statistics (m x n matrix)
   */
SEXP GetODPResult(SEXP tid)
{
  int nProt = 0;
  PROTECT(tid = AS_INTEGER(tid)); nProt++;
  THREAD_ID_T ttid = (THREAD_ID_T)INTEGER_VALUE(tid);

  COdpThread *pThread = (COdpThread*)CProgressThreadManager::GetProgressThread(ttid);

  if(NULL == pThread)
    {
      Rf_error("Could not get thread information for thread id = %d", ttid);
      UNPROTECT(nProt);
      return R_NilValue;
    }

  if(pThread->GetIsCanceled())
    {
      Rf_error("Thread was canceled");
      UNPROTECT(nProt);
      return R_NilValue;
    }
  if(!pThread->GetIsFinished())
    {
      Rf_error("Thread is still running");
      UNPROTECT(nProt);
      return R_NilValue;
    }
  
  /* Create a list with "lr" and "lr0" members */
  SEXP slist;
  PROTECT(slist = GetODPResultList(pThread));nProt++;
  UNPROTECT(nProt);
  return(slist);
}

/*
 * Constructs an R-style list containing
 * ODP results.
 */

SEXP GetODPResultList (COdpThread *pThread)
{
  int nProt = 0;

  if (NULL == pThread)
    {
      Rf_error("null pThread passed to GetODPResultList");
      return R_NilValue;
    }

  SEXP slist, snames, sval;
  int i = 0;
  PROTECT(slist = allocVector(VECSXP, 2)); nProt++;
  PROTECT(snames = allocVector(STRSXP, 2)); nProt++;
  
  int size = pThread->GetObservedResultSize();
  double* pObservedData = pThread->GetObservedResult();
  PROTECT(sval=NEW_NUMERIC(size)); nProt++;
  double* pResult = NUMERIC_POINTER(sval);
  for(int j = 0; j < size ; j++)
    {
      pResult[j] = pObservedData[j];
    }
  SET_STRING_ELT(snames, i, CREATE_STRING_VECTOR("lr"));
  SET_VECTOR_ELT(slist, i, sval);
  i++;

  size = pThread->GetNullResultSize();
  double* pNullData = pThread->GetNullResult();
  PROTECT(sval=NEW_NUMERIC(size)); nProt++;
  pResult = NUMERIC_POINTER(sval);
  for(int j = 0; j < size; j++)
    {
      pResult[j] = pNullData[j];
    }
  SET_STRING_ELT(snames, i, CREATE_STRING_VECTOR("lr0"));
  SET_VECTOR_ELT(slist, i, sval);
  i++;

  setAttrib(slist, R_NamesSymbol, snames);

  UNPROTECT(nProt);
  return slist;

}

  /*
   * Gets the current version of the EDGE library.
   */
SEXP GetEdgeLibraryVersion()
{
  SEXP sversion;

  PROTECT(sversion = NEW_STRING(1));
  SET_STRING_ELT(sversion, 0, mkChar(EDGE_LIBRARY_VERSION));
  UNPROTECT(1);
  return(sversion);
}

  /* 
   * Starts the "intercept" version of the timecourse function  
   * in a background thread
   * See corresponding function in edge.c for clues about what the arguments are
   */
// SEXP RunTimecourseInt(SEXP fit, SEXP res, SEXP v, 
// 			SEXP v1, SEXP H1, SEXP H0, 
// 			SEXP alpha1, SEXP ind, SEXP idx, 
// 			SEXP eps, SEXP m, SEXP n, SEXP N, 
// 			SEXP nGrp, SEXP B,
// 			// threading?
// 			SEXP startThread
// 			)
SEXP RunTimecourseInt(SEXP arglist, SEXP startThread)
{
  int nProtected = 0;
  
  SEXP tid;
  PROTECT(tid = NEW_INTEGER(1));  nProtected++;

//  PROTECT(fit = AS_NUMERIC(fit));  nProtected++;
//  PROTECT(res = AS_NUMERIC(res));  nProtected++;
//  PROTECT(v = AS_INTEGER(v));  nProtected++;
//  PROTECT(v1 = AS_INTEGER(v1));  nProtected++;
//  PROTECT(H1 = AS_NUMERIC(H1));  nProtected++;
//  PROTECT(H0 = AS_NUMERIC(H0));  nProtected++;
//  PROTECT(alpha1 = AS_NUMERIC(alpha1));  nProtected++;
//  PROTECT(ind = AS_INTEGER(ind));  nProtected++;
//  PROTECT(idx = AS_INTEGER(idx));  nProtected++;
//  PROTECT(eps = AS_NUMERIC(eps));  nProtected++;
//  PROTECT(m = AS_INTEGER(m));  nProtected++;
//  PROTECT(n = AS_INTEGER(n));  nProtected++;
//  PROTECT(N = AS_INTEGER(N));  nProtected++;
//  PROTECT(nGrp = AS_INTEGER(nGrp));  nProtected++;
//  PROTECT(B = AS_INTEGER(B));  nProtected++;
//  PROTECT(startThread = AS_LOGICAL(startThread)); nProtected++;
//
//  double *xfit = NUMERIC_POINTER(fit);
//  double *xres = NUMERIC_POINTER(res);
//  int *xv = INTEGER_POINTER(v);
//  int *xv1 = INTEGER_POINTER(v1);
//  double *xH1 = NUMERIC_POINTER(H1);
//  double *xH0 = NUMERIC_POINTER(H0);
//  double *xalpha1 = NUMERIC_POINTER(alpha1);
//  int *xind = INTEGER_POINTER(ind);
//  int *xidx = INTEGER_POINTER(idx);
//  double xeps = NUMERIC_VALUE(eps);
//  long xm = INTEGER_VALUE(m);
//  long xn = INTEGER_VALUE(n);
//  long xN = INTEGER_VALUE(N);
// int *xnGrp = INTEGER_POINTER(nGrp);
// long xB = INTEGER_VALUE(B);
  bool bStartThread = LOGICAL_POINTER(startThread)[0];

  /* Sanity checks */
// move checks to constructor?
//  CheckLength(GET_LENGTH(fit), "fit", xm*xn);
//  CheckLength(GET_LENGTH(res), "res", xm*xn);
//  CheckLength(GET_LENGTH(v), "v", xn*xB);
//  CheckLength(GET_LENGTH(v1), "v1", xn*xB);
//  CheckLength(GET_LENGTH(H1), "H1", xn*xn);
//  CheckLength(GET_LENGTH(H0), "H0", xn*xn);
//  CheckLength(GET_LENGTH(alpha1), "alpha1", xm*xn);
//  CheckLength(GET_LENGTH(ind), "ind", xn);
//  CheckLength(GET_LENGTH(idx), "idx", xn);
//  CheckLength(GET_LENGTH(nGrp), "nGrp", xN);

  CTimecourseIntercept* pThread = NULL;
  try {
    pThread = new CTimecourseIntercept(arglist);
  }
  catch(bad_alloc&) {
    pThread = NULL;
  }

  if(NULL == pThread)
    {
      Rf_error("Could not allocate CTimecourseThread object");
      UNPROTECT(nProtected);
      return R_NilValue;
    }

  int status;
  if (!bStartThread)
    {
      /* Do computation and return when finished */
      status = pThread->Go(0);
      if(0 != status)
	{
	  /* Computation failed; cleanup and return */
	  delete pThread;
	  Rf_error("Timecourse computation failed (error number %d)", status);
	  UNPROTECT(nProtected);
	  return R_NilValue;
	}
      else
	{
	  /* Get result data and return */
	  SEXP slist;
	  PROTECT(slist = GetTimecourseResultList(pThread)); nProtected++;
	  delete pThread;
	  UNPROTECT(nProtected);
	  return(slist);
	}
    }
  else
    {
      status = pThread->Start();
      if(0 != status)
	{
	  if(!pThread->GetIsStarted())
	    {
	      delete pThread;
	      Rf_error("Could not start thread (error number %d)", status);
	      UNPROTECT(nProtected);
	      return R_NilValue;
	    }
	  else
	    Rf_warning("Thread %d was started, but returned status code %d", 
		       pThread->GetThreadId(), status);
	}
      
      THREAD_ID_T xtid = pThread->GetThreadId();
      INTEGER_POINTER(tid)[0] = (int)xtid;
      
      UNPROTECT(nProtected);
      return tid;
    }
}

  /*
   * Starts a background thread that calls "nullFLongNoInt".
   * Arguments:
   *  fit0: Vector of length nrow(H0)*nrow(dat)
   *  res: Vector of length m*n
   *  v: Vector length n*B
   *  rmv: vector length N
   *  inc: n-N vector
   *  H1: Vector same length as H0
   *  H0: Vector length m*(n/2)
   *  wsq: vector length n (or is it n-N * n-N)
   *  ind: vector length n
   *  idx: vector length n
   *  B: Integer
   *  m: number of rows in dat
   *  n: number of columns in dat
   *  N: integer
   *  nGrp: N vector
   *  
   */
// SEXP RunTimecourseNoInt(SEXP fit, SEXP res, SEXP v, 
// 			  SEXP rmv, SEXP inc, SEXP H1, SEXP H0,
// 			  SEXP wsq, SEXP ind, SEXP idx, SEXP m, 
// 			  SEXP n, SEXP N, SEXP nGrp, SEXP B,
// 			  // threading?
// 			  SEXP startThread
// 			  )
SEXP RunTimecourseNoInt(SEXP arglist, SEXP startThread)
{
  int nProtected = 0;
  
  SEXP tid;
  PROTECT(tid = NEW_INTEGER(1));  nProtected++;

//  PROTECT(fit = AS_NUMERIC(fit));  nProtected++;
//  PROTECT(res = AS_NUMERIC(res));  nProtected++;
//  PROTECT(v = AS_INTEGER(v));  nProtected++;
//  PROTECT(rmv = AS_INTEGER(rmv));  nProtected++;
//  PROTECT(inc = AS_INTEGER(inc));  nProtected++;
//  PROTECT(H1 = AS_NUMERIC(H1));  nProtected++;
//  PROTECT(H0 = AS_NUMERIC(H0));  nProtected++;
//  PROTECT(wsq = AS_NUMERIC(wsq));  nProtected++;
//  PROTECT(ind = AS_INTEGER(ind));  nProtected++;
//  PROTECT(idx = AS_INTEGER(idx));  nProtected++;
//  PROTECT(m = AS_INTEGER(m));  nProtected++;
//  PROTECT(n = AS_INTEGER(n));  nProtected++;
//  PROTECT(N = AS_INTEGER(N));  nProtected++;
//  PROTECT(nGrp = AS_INTEGER(nGrp));  nProtected++;
//  PROTECT(B = AS_INTEGER(B));  nProtected++;
//  PROTECT(startThread = AS_LOGICAL(startThread)); nProtected++;

//  double *xfit = NUMERIC_POINTER(fit);
//  double *xres = NUMERIC_POINTER(res);
//  int *xv = INTEGER_POINTER(v);
//  int *xrmv = INTEGER_POINTER(rmv);
//  int *xinc = INTEGER_POINTER(inc);
//  double *xH1 = NUMERIC_POINTER(H1);
//  double *xH0 = NUMERIC_POINTER(H0);
//  double *xwsq = NUMERIC_POINTER(wsq);
//  int *xind = INTEGER_POINTER(ind);
//  int *xidx = INTEGER_POINTER(idx);
//  long xm = INTEGER_VALUE(m);
//  long xn = INTEGER_VALUE(n);
//  long xN = INTEGER_VALUE(N);
//  int *xnGrp = INTEGER_POINTER(nGrp);
//  long xB = INTEGER_VALUE(B);
  bool bStartThread = LOGICAL_POINTER(startThread)[0];

  /* Sanity checks */
// gah - Sanity checks moved to TimecourseNoInt.cpp

//  CheckLength(GET_LENGTH(fit), "fit", xm*xn);
//  CheckLength(GET_LENGTH(res), "res", xm*(xn - xN));
//  CheckLength(GET_LENGTH(v), "v", (xn-xN)*xB);
//  CheckLength(GET_LENGTH(rmv), "rmv", xN);
//  CheckLength(GET_LENGTH(inc), "inc", xn-xN);
//  CheckLength(GET_LENGTH(H1), "H1", xn*xn);
//  CheckLength(GET_LENGTH(H0), "H0", xn*xn);
//  CheckLength(GET_LENGTH(wsq), "wsq", (xn-xN)*(xn-xN));
//  CheckLength(GET_LENGTH(ind), "ind", xn);
//  CheckLength(GET_LENGTH(idx), "idx", xn);
//  CheckLength(GET_LENGTH(nGrp), "nGrp", xN);

  CTimecourseNoIntercept* pThread = NULL;
  try {
//    pThread = new CTimecourseNoIntercept(xfit, xres, xv, xrmv, xinc, 
//			       xH1, xH0, xwsq, xind, xidx, 
//			       xm, xn, xN, xnGrp, xB);
    pThread = new CTimecourseNoIntercept(arglist);
  }
  catch(bad_alloc&) {
    pThread = NULL;
  }

  if(NULL == pThread)
    {
      Rf_error("Could not allocate CTimecourseThread object");
      UNPROTECT(nProtected);
      return R_NilValue;
    }

  int status;
  if (!bStartThread)
    {
      /* Do computation and return when finished */
      status = pThread->Go(0);
      if(0 != status)
	{
	  /* Computation failed; cleanup and return */
	  delete pThread;
	  Rf_error("Timecourse computation failed (error number %d)", status);
	  UNPROTECT(nProtected);
	  return R_NilValue;
	}
      else
	{
	  /* Get result data and return */
	  SEXP slist;
	  PROTECT(slist = GetTimecourseResultList(pThread)); nProtected++;
	  delete pThread;
	  UNPROTECT(nProtected);
	  return(slist);
	}
    }
  else
    {
      status = pThread->Start();
      if(0 != status)
	{
	  if(!pThread->GetIsStarted())
	    {
	      delete pThread;
	      Rf_error("Could not start thread (error number %d)", status);
	      UNPROTECT(nProtected);
	      return R_NilValue;
	    }
	  else
	    Rf_warning("Thread %d was started, but returned status code %d", 
		       pThread->GetThreadId(), status);
	}
      
      THREAD_ID_T xtid = pThread->GetThreadId();
      INTEGER_POINTER(tid)[0] = (int)xtid;
      
      UNPROTECT(nProtected);
      return tid;
    }
}

  /*
   * Starts a background thread that calls "nullFIndep".
   * Arguments:
   *  fit0: Vector of length nrow(H0)*nrow(dat)
   *  res: Vector of length nrow(dat)*ncol(dat)
   *  v: Vector length n*B
   *  H0: Vector length nrow(dat)*ncol(dat)/2
   *  H1: Vector same length as H0
   *  B: Integer
   *  m: number of rows in dat
   *  n: number of columns in dat
   */
// SEXP RunTimecourseIndep(SEXP fit0, SEXP res, SEXP v, SEXP H1, 
// 			  SEXP H0, SEXP B, SEXP m, SEXP n,
// 			  // threading?
// 			  SEXP startThread
// 			  )
SEXP RunTimecourseIndep(SEXP arglist, SEXP startThread)
{
  int nProtected = 0;
  
  SEXP tid;
  PROTECT(tid = NEW_INTEGER(1));  nProtected++;

//  PROTECT(fit0 = AS_NUMERIC(fit0));  nProtected++;
//  PROTECT(res = AS_NUMERIC(res));  nProtected++;
//  PROTECT(v = AS_INTEGER(v));  nProtected++;
//  PROTECT(H1 = AS_NUMERIC(H1));  nProtected++;
//  PROTECT(H0 = AS_NUMERIC(H0));  nProtected++;
//  PROTECT(B = AS_INTEGER(B));  nProtected++;
//  PROTECT(m = AS_INTEGER(m));  nProtected++;
//  PROTECT(n = AS_INTEGER(n));  nProtected++;
//  PROTECT(startThread = AS_LOGICAL(startThread)); nProtected++;

//  double *xfit0 = NUMERIC_POINTER(fit0);
//  double *xres = NUMERIC_POINTER(res);
//  int *xv = INTEGER_POINTER(v);
//  double *xH1 = NUMERIC_POINTER(H1);
//  double *xH0 = NUMERIC_POINTER(H0);
//  int xB = INTEGER_VALUE(B);
//  int xm = INTEGER_VALUE(m);
//  int xn = INTEGER_VALUE(n);
  bool bStartThread = LOGICAL_POINTER(startThread)[0];

  /* Sanity checks */
//  CheckLength(GET_LENGTH(fit0), "fit0", xm*xn);
//  CheckLength(GET_LENGTH(res), "res", xm*xn);
//  CheckLength(GET_LENGTH(v), "v", xn*xB);
//  CheckLength(GET_LENGTH(H0), "H0", xn*xn);
//  CheckLength(GET_LENGTH(H1), "H1", GET_LENGTH(H0));
//  
  CTimecourseThread* pThread = NULL;
  try { 
//    pThread =new CTimecourseThread(xfit0, xres, xv, xH1, xH0, xB, xm, xn);
    pThread =new CTimecourseThread(arglist);
  }
  catch(bad_alloc&) {
    pThread = NULL;
  }

  if(NULL == pThread)
    {
      Rf_error("Could not allocate CTimecourseThread object");
      UNPROTECT(nProtected);
      return R_NilValue;
    }

  int status;
  if (!bStartThread)
    {
      /* Do computation and return when finished */
      status = pThread->Go(0);
      if(0 != status)
	{
	  /* Computation failed; cleanup and return */
	  delete pThread;
	  Rf_error("Timecourse computation failed (error number %d)", status);
	  UNPROTECT(nProtected);
	  return R_NilValue;
	}
      else
	{
	  /* Get result data and return */
	  SEXP slist;
	  PROTECT(slist = GetTimecourseResultList(pThread)); nProtected++;
	  delete pThread;
	  UNPROTECT(nProtected);
	  return(slist);
	}
    }
  else
    {
      status = pThread->Start();
      if(0 != status)
	{
	  if(!pThread->GetIsStarted())
	    {
	      delete pThread;
	      Rf_error("Could not start thread (error number %d)", status);
	      UNPROTECT(nProtected);
	      return R_NilValue;
	    }
	  else
	    Rf_warning("Thread %d was started, but returned status code %d", 
		       pThread->GetThreadId(), status);
	}
      
      THREAD_ID_T xtid = pThread->GetThreadId();
      INTEGER_POINTER(tid)[0] = (int)xtid;
      
      UNPROTECT(nProtected);
      return tid;
    }
}

SEXP GetTimecourseResult(SEXP tid)
{
  int nProt = 0;
  PROTECT(tid = AS_INTEGER(tid)); nProt++;
  THREAD_ID_T ttid = (THREAD_ID_T)INTEGER_VALUE(tid);

  CProgressThread *pThread = (CProgressThread*)
    CProgressThreadManager::GetProgressThread(ttid);

  if(NULL == pThread)
    {
      Rf_error("Could not get thread information for thread id = %d", ttid);
      UNPROTECT(nProt);
      return R_NilValue;
    }

  if(pThread->GetIsCanceled())
    {
      Rf_error("Thread was canceled");
      return R_NilValue;
    }
  if(!pThread->GetIsFinished())
    {
      Rf_error("Thread is still running");
      return R_NilValue;
    }
  
  SEXP slist;
  PROTECT(slist = GetTimecourseResultList(pThread)); nProt++;
  UNPROTECT(nProt);

  return slist;
}

SEXP GetTimecourseResultList(CProgressThread *pThread)
{
  int nProt = 0;

  if (NULL == pThread)
    {
      Rf_error("Null argument passed to GetTimecourseResultList");
      return R_NilValue;
    }

  // Create a list with results
  SEXP slist, snames, sval;
  int i = 0;
  PROTECT(slist = allocVector(VECSXP, 2)); nProt++;
  PROTECT(snames = allocVector(STRSXP, 2)); nProt++;
  
  int size = pThread->GetReturnDataSize();
  if(size < 0)
    {
      Rf_error("Invalid return data size: %d", size);
      UNPROTECT(nProt);
      return R_NilValue;
    }

  double* pResultArray = (double*)pThread->GetReturnData();
  if(NULL == pResultArray)
    {
      Rf_error("Result data not found");
      UNPROTECT(nProt);
      return R_NilValue;
    }

  PROTECT(sval=NEW_NUMERIC(size)); nProt++;
  double* pResultSexp = NUMERIC_POINTER(sval);
  for(int j = 0; j < size ; j++)
    {
      pResultSexp[j] = pResultArray[j];
    }
  SET_STRING_ELT(snames, i, CREATE_STRING_VECTOR("F"));
  SET_VECTOR_ELT(slist, i, sval);
  i++;

  // Associate name strings with list values
  setAttrib(slist, R_NamesSymbol, snames);

  UNPROTECT(nProt);
  return slist;

}

#ifdef __cplusplus
}
#endif
