C++ INTEGER函数代码示例

void getDrivesForArray(const std::string& arrayName,
                       MDInterface& md,
                       QueryData& data) {
  std::string path(md.getPathByDevName(arrayName));
  if (path.empty()) {
    LOG(ERROR) << "Could not get file path for " << arrayName;
    return;
  }

  mdu_array_info_t array;
  if (!md.getArrayInfo(path, array)) {
    return;
  }

  /* Create a vector of with all expected slot positions.  As we work through
   * the RAID disks, we remove discovered slots */
  std::vector<size_t> missingSlots(array.raid_disks);
  std::iota(missingSlots.begin(), missingSlots.end(), 0);

  /* Keep track of index in QueryData that have removed slots since we can't
   * make safe assumptions about it's original slot position if disk_number >=
   * total_disk and we're unable to deteremine total number of missing slots
   * until we walk thru all MD_SB_DISKS */
  std::vector<size_t> removedSlots;

  size_t qdPos = data.size();
  for (size_t i = 0; i < MD_SB_DISKS; i++) {
    mdu_disk_info_t disk;
    disk.number = i;
    if (!md.getDiskInfo(path, disk)) {
      continue;
    }

    if (disk.major > 0) {
      Row r;
      r["md_device_name"] = arrayName;
      r["drive_name"] = md.getDevName(disk.major, disk.minor);
      r["state"] = getDiskStateStr(disk.state);

      if (disk.raid_disk >= 0) {
        r["slot"] = INTEGER(disk.raid_disk);
        missingSlots.erase(
            std::remove(
                missingSlots.begin(), missingSlots.end(), disk.raid_disk),
            missingSlots.end());

        /* We assume that if the disk number is less than the total disk count
         * of the array, then it assumes its original slot position;  If the
         * number is greater than the disk count, then it's not safe to make
         * that assumption. We do this check here b/c if a recovery is targeted
         * for the same slot, we potentially miss identifying the original slot
         * position of the bad disk. */
      } else if (disk.raid_disk < 0 && disk.number < array.raid_disks) {
        r["slot"] = std::to_string(disk.number);
        missingSlots.erase(
            std::remove(missingSlots.begin(), missingSlots.end(), disk.number),
            missingSlots.end());

        /* Mark QueryData position as a removedSlot to handle later*/
      } else {
        removedSlots.push_back(qdPos);
      }

      qdPos++;
      data.push_back(r);
    }
  }

  /* Handle all missing slots.  See `scattered_faulty_and_removed` unit test in
   * `./tests/md_tables_tests.cpp`*/
  for (const auto& slot : missingSlots) {
    if (!removedSlots.empty()) {
      data[removedSlots[0]]["slot"] = INTEGER(slot);
      removedSlots.erase(removedSlots.begin());

    } else {
      Row r;
      r["md_device_name"] = arrayName;
      r["drive_name"] = "unknown";
      r["state"] = "removed";
      r["slot"] = std::to_string(slot);
      data.push_back(r);
    }
  }
}

int * AllocProtectInt(int n){

  SEXP x = R_NilValue;
  PROTECT(x = NEW_INTEGER(n));
  return( INTEGER(x) );
}

SEXP agmart3(SEXP surv2, SEXP score2, SEXP weight2, SEXP strata2,
	     SEXP sortx, SEXP method2) {

    int k, ksave;
    int p, istrat, indx2;
    double deaths, denom, e_denom;
    double hazard, e_hazard, cumhaz;
    double temp, time;
    double wtsum;
    int n, person;
 
    /* pointers to the input data */
    double *start, *stop, *event;
    double *weight, *score;
    int    *sort1, *sort2, *strata;

    int method;  /* integer version of input */

    /* output */
    SEXP resid2;
    double *resid;

    n = nrows(surv2);
    method = asInteger(method2);
    start = REAL(surv2);
    stop  = start +n;
    event = stop +n;
    weight= REAL(weight2);
    score = REAL(score2);
    sort1 = INTEGER(sortx);
    sort2 = sort1 + n;
    strata= INTEGER(strata2);

    PROTECT(resid2 = allocVector(REALSXP, n));
    resid = REAL(resid2);

    /*
    **  'person' walks through the the data from 1 to n,
    **     sort1[0] points to the largest stop time, sort1[1] the next, ...
    **  'time' is a scratch variable holding the time of current interest
    **  'indx2' walks through the start times.  It will be smaller than 
    **    'person': if person=27 that means that 27 subjects have stop >=time,
    **    and are thus potential members of the risk set.  If 'indx2' =9,
    **    that means that 9 subjects have start >=time and thus are NOT part
    **    of the risk set.  (stop > start for each subject guarrantees that
    **    the 9 are a subset of the 27). 
    **  Basic algorithm: move 'person' forward, adding the new subject into
    **    the risk set.  If this is a new, unique death time, take selected
    **    old obs out of the sums, add in obs tied at this time, then update
    **    the cumulative hazard. Everything resets at the end of a stratum.
    **  The sort order is from large time to small, so we encounter a subject's
    **    ending time first, then their start time.
    **  The martingale residual for a subject is 
    **     status - (cumhaz at end of their interval - cumhaz at start)*score
    */
    istrat=0;
    indx2 =0;
    denom =0;
    cumhaz =0;
    for (person=0; person <n; ) {
	p = sort1[person];
	if (event[p] ==0) { /* censored */
	    denom += score[p] * weight[p];
	    resid[p] = cumhaz * score[p];
	    person++;
	} else {
	    time = stop[p]; /* found a new, unique death time */
	    /* 
	    ** Remove those subjects whose start time is to the right
	    **  from the risk set, and finish computation of their residual
	    */
	    for (;  indx2 <strata[istrat]; indx2++) {
		p = sort2[indx2];
		if (start[p] < time) break;
		denom -= score[p] * weight[p];
		resid[p] -= cumhaz * score[p];
	    }

	    /*
	    **	Add up over this death time, for all subjects
	    */
	    deaths =0;
	    e_denom =0;
	    wtsum =0;
	    for (k=person; k<strata[istrat]; k++) {
		p = sort1[k];
		if (stop[p]  < time) break;  /* only tied times */ 
		denom += score[p] * weight[p];
		if (event[p] ==1) {
		    deaths ++;
		    e_denom += score[p] * weight[p];
		    wtsum += weight[p];
		}
	    }
	    ksave = k;
	    
	    /* compute the increment in hazard 
	    ** hazard = usual increment
	    ** e_hazard = efron increment, for tied deaths only
	    */
//.........这里部分代码省略.........

SEXP gbm_pred
(
   SEXP radX,        // the data matrix
   SEXP rcRows,      // number of rows
   SEXP rcCols,      // number of columns
   SEXP rcTrees,     // number of trees, may be a vector
   SEXP rdInitF,     // the initial value
   SEXP rTrees,      // the list of trees
   SEXP rCSplits,    // the list of categorical splits
   SEXP raiVarType,  // indicator of continuous/nominal
   SEXP riSingleTree // boolean whether to return only results for one tree
)
{
   unsigned long hr = 0;
   int iTree = 0;
   int iObs = 0;
   int cRows = INTEGER(rcRows)[0];
   int cPredIterations = LENGTH(rcTrees);
   int iPredIteration = 0;
   int cTrees = 0;

   SEXP rThisTree = NULL;
   int *aiSplitVar = NULL;
   double *adSplitCode = NULL;
   int *aiLeftNode = NULL;
   int *aiRightNode = NULL;
   int *aiMissingNode = NULL;
   int iCurrentNode = 0;
   double dX = 0.0;
   int iCatSplitIndicator = 0;
   bool fSingleTree = (INTEGER(riSingleTree)[0]==1);

   SEXP radPredF = NULL;

   // allocate the predictions to return
   PROTECT(radPredF = allocVector(REALSXP, cRows*cPredIterations));
   if(radPredF == NULL)
   {
      hr = GBM_OUTOFMEMORY;
      goto Error;
   }

   // initialize the predicted values
   if(!fSingleTree)
   {
      // initialize with the intercept for only the smallest rcTrees
      for(iObs=0; iObs<cRows; iObs++)
      {
         REAL(radPredF)[iObs] = REAL(rdInitF)[0];
      }
   }
   else
   {
      for(iObs=0; iObs<cRows*cPredIterations; iObs++)
      {
         REAL(radPredF)[iObs] = 0.0;
      }
   }

   iTree = 0;
   for(iPredIteration=0; iPredIteration<LENGTH(rcTrees); iPredIteration++)
   {
      cTrees = INTEGER(rcTrees)[iPredIteration];
      if(fSingleTree) iTree=cTrees-1;
      if(!fSingleTree && (iPredIteration>0))
      {
         // copy over from the last rcTrees
         for(iObs=0; iObs<cRows; iObs++)
         {
            REAL(radPredF)[cRows*iPredIteration+iObs] =
               REAL(radPredF)[cRows*(iPredIteration-1)+iObs];
         }
      }
      while(iTree<cTrees)
      {
         rThisTree     = VECTOR_ELT(rTrees,iTree);
         // these relate to columns returned by pretty.gbm.tree()
         aiSplitVar    = INTEGER(VECTOR_ELT(rThisTree,0));
         adSplitCode   = REAL   (VECTOR_ELT(rThisTree,1));
         aiLeftNode    = INTEGER(VECTOR_ELT(rThisTree,2));
         aiRightNode   = INTEGER(VECTOR_ELT(rThisTree,3));
         aiMissingNode = INTEGER(VECTOR_ELT(rThisTree,4));
         for(iObs=0; iObs<cRows; iObs++)
         {
            iCurrentNode = 0;
            while(aiSplitVar[iCurrentNode] != -1)
            {
               dX = REAL(radX)[aiSplitVar[iCurrentNode]*cRows + iObs];
               // missing?
               if(ISNA(dX))
               {
                  iCurrentNode = aiMissingNode[iCurrentNode];
               }
               // continuous?
               else if(INTEGER(raiVarType)[aiSplitVar[iCurrentNode]] == 0)
               {
                  if(dX < adSplitCode[iCurrentNode])
                  {
                        iCurrentNode = aiLeftNode[iCurrentNode];
                  }
//.........这里部分代码省略.........

SEXP gbm
(
    SEXP radY,       // outcome or response
    SEXP radOffset,  // offset for f(x), NA for no offset
    SEXP radX,        
    SEXP raiXOrder,        
    SEXP radWeight,
    SEXP radMisc,   // other row specific data (eg failure time), NA=no Misc
    SEXP rcRows,
    SEXP rcCols,
    SEXP racVarClasses,
    SEXP ralMonotoneVar,
    SEXP rszFamily, 
    SEXP rcTrees,
    SEXP rcDepth,       // interaction depth
    SEXP rcMinObsInNode,
    SEXP rdShrinkage,
    SEXP rdBagFraction,
    SEXP rcTrain,
    SEXP radFOld,
    SEXP rcCatSplitsOld,
    SEXP rcTreesOld,
    SEXP rfVerbose
)
{
    unsigned long hr = 0;

    SEXP rAns = NULL;
    SEXP rNewTree = NULL;
    SEXP riSplitVar = NULL;
    SEXP rdSplitPoint = NULL;
    SEXP riLeftNode = NULL;
    SEXP riRightNode = NULL;
    SEXP riMissingNode = NULL;
    SEXP rdErrorReduction = NULL;
    SEXP rdWeight = NULL;
    SEXP rdPred = NULL;

    SEXP rdInitF = NULL;
    SEXP radF = NULL;
    SEXP radTrainError = NULL;
    SEXP radValidError = NULL;
    SEXP radOOBagImprove = NULL;

    SEXP rSetOfTrees = NULL;
    SEXP rSetSplitCodes = NULL;
    SEXP rSplitCode = NULL;

    VEC_VEC_CATEGORIES vecSplitCodes;

    int i = 0;
    int iT = 0;
    int cTrees = INTEGER(rcTrees)[0];
    const int cResultComponents = 7;
    // rdInitF, radF, radTrainError, radValidError, radOOBagImprove
    // rSetOfTrees, rSetSplitCodes
    const int cTreeComponents = 8;
    // riSplitVar, rdSplitPoint, riLeftNode,
    // riRightNode, riMissingNode, rdErrorReduction, rdWeight, rdPred
    int cNodes = 0;
    int cTrain = INTEGER(rcTrain)[0];

    double dTrainError = 0.0;
    double dValidError = 0.0;
    double dOOBagImprove = 0.0;

    CGBM *pGBM = NULL;
    CDataset *pData = NULL;
    CDistribution *pDist = NULL;

    // set up the dataset
    pData = new CDataset();
    if(pData==NULL)
    {
        hr = GBM_OUTOFMEMORY;
        goto Error;
    }

    // initialize R's random number generator
    GetRNGstate();

    // initialize some things
    hr = gbm_setup(REAL(radY),
                   REAL(radOffset),
                   REAL(radX),
                   INTEGER(raiXOrder),
                   REAL(radWeight),
                   REAL(radMisc),
                   INTEGER(rcRows)[0],
                   INTEGER(rcCols)[0],
                   INTEGER(racVarClasses),
                   INTEGER(ralMonotoneVar),
                   CHAR(STRING_ELT(rszFamily,0)),
                   INTEGER(rcTrees)[0],
                   INTEGER(rcDepth)[0],
                   INTEGER(rcMinObsInNode)[0],
                   REAL(rdShrinkage)[0],
                   REAL(rdBagFraction)[0],
                   INTEGER(rcTrain)[0],
                   pData,
//.........这里部分代码省略.........

  SEXP MBAPoints(SEXP xyz, SEXP xyzEst, SEXP m, SEXP n, SEXP h, SEXP extend, SEXP verbose) {

    int i,j;
    int nProtect = 0;
    SEXP xyzPts, Z;

    //get the surface points
    PROTECT(xyzPts = getAttrib(xyz, R_DimSymbol)); nProtect++;      
    int nPts = INTEGER(xyzPts)[0];

    typedef std::vector<double> dVec;
    boost::shared_ptr<dVec> x_arr(new std::vector<double>);
    boost::shared_ptr<dVec> y_arr(new std::vector<double>);
    boost::shared_ptr<dVec> z_arr(new std::vector<double>);

    for(i = 0; i < nPts; i++){
      x_arr->push_back(REAL(xyz)[i]);
      y_arr->push_back(REAL(xyz)[nPts+i]);
      z_arr->push_back(REAL(xyz)[2*nPts+i]);
    }

    double maxXSurf = *std::max_element((*x_arr).begin(), (*x_arr).end());
    double minXSurf = *std::min_element((*x_arr).begin(), (*x_arr).end());
    double maxYSurf = *std::max_element((*y_arr).begin(), (*y_arr).end());
    double minYSurf = *std::min_element((*y_arr).begin(), (*y_arr).end());

    //get the points
    xyzPts = getAttrib(xyzEst, R_DimSymbol);      
    int nEstPts = INTEGER(xyzPts)[0];

    double maxXPt = REAL(xyzEst)[0];
    double minXPt = REAL(xyzEst)[0];
    double maxYPt = REAL(xyzEst)[nEstPts];
    double minYPt = REAL(xyzEst)[nEstPts];

    for(i = 0; i < nEstPts; i++){
      if(REAL(xyzEst)[i] > maxXPt) maxXPt = REAL(xyzEst)[i];
      if(REAL(xyzEst)[i] < minXPt) minXPt = REAL(xyzEst)[i];
      if(REAL(xyzEst)[nEstPts+i] > maxYPt) maxYPt = REAL(xyzEst)[nEstPts+i];
      if(REAL(xyzEst)[nEstPts+i] < minYPt) minYPt = REAL(xyzEst)[nEstPts+i];
    }

    std::string extendDirection = "";
    if(INTEGER(extend)[0]){
      if(maxXSurf < maxXPt){
	maxXSurf = maxXPt;
	extendDirection="+x ";
      }
      if(minXSurf > minXPt){
	minXSurf = minXPt;
	extendDirection+="-x ";
      }
      if(maxYSurf < maxYPt){
	maxYSurf = maxYPt;
	extendDirection+="+y ";
      }
      if(minYSurf > minYPt){ 
	minYSurf = minYPt;
	extendDirection+="-y ";
      }
    }

    //init
    MBA mba(x_arr, y_arr, z_arr);
    
    mba.setDomain(minXSurf, minYSurf, maxXSurf, maxYSurf);

    mba.MBAalg(INTEGER(m)[0], INTEGER(n)[0], INTEGER(h)[0]);

    //retrieve the spline surface and evaluate
    UCBspl::SplineSurface surface = mba.getSplineSurface(); 

    double umin, vmin, umax, vmax;
    surface.getDomain(umin, vmin, umax, vmax);
    
    PROTECT(Z = allocVector(REALSXP, nEstPts)); nProtect++;
    
    int ptsOutSide = 0;
    for (i = 0; i < nEstPts; i++){
      if(REAL(xyzEst)[i] < umin || REAL(xyzEst)[i] > umax ||
	 REAL(xyzEst)[nEstPts+i] < vmin || REAL(xyzEst)[nEstPts+i] > vmax){
	REAL(Z)[i] = NA_REAL;
	ptsOutSide++;
      }else{ 
	REAL(Z)[i] = surface.f(REAL(xyzEst)[i], REAL(xyzEst)[nEstPts+i]);
      }
    }

    if(INTEGER(verbose)[0]){
      if(extendDirection != "")
	warning("domain extended in the %sdirection(s)\n", extendDirection.c_str());
      
      if(ptsOutSide)
	warning("%i point(s) fell outside the domain and were set to NA\n", ptsOutSide);
    }

    //clean-up
    mba.cleanup(2);

    //just to be sure
//.........这里部分代码省略.........

static void extractItem(char *buffer, SEXP ans, int i, LocalData *d)
{
    char *endp;
    switch(TYPEOF(ans)) {
    case NILSXP:
	break;
    case LGLSXP:
	if (isNAstring(buffer, 0, d))
	    LOGICAL(ans)[i] = NA_INTEGER;
	else {
	    int tr = StringTrue(buffer), fa = StringFalse(buffer);
	    if(tr || fa) LOGICAL(ans)[i] = tr;
	    else expected("a logical", buffer, d);
	}
	break;
    case INTSXP:
	if (isNAstring(buffer, 0, d))
	    INTEGER(ans)[i] = NA_INTEGER;
	else {
	    INTEGER(ans)[i] = Strtoi(buffer, 10);
	    if (INTEGER(ans)[i] == NA_INTEGER)
		expected("an integer", buffer, d);
	}
	break;
    case REALSXP:
	if (isNAstring(buffer, 0, d))
	    REAL(ans)[i] = NA_REAL;
	else {
	    REAL(ans)[i] = Strtod(buffer, &endp, TRUE, d);
	    if (!isBlankString(endp))
		expected("a real", buffer, d);
	}
	break;
    case CPLXSXP:
	if (isNAstring(buffer, 0, d))
	    COMPLEX(ans)[i].r = COMPLEX(ans)[i].i = NA_REAL;
	else {
	    COMPLEX(ans)[i] = strtoc(buffer, &endp, TRUE, d);
	    if (!isBlankString(endp))
		expected("a complex", buffer, d);
	}
	break;
    case STRSXP:
	if (isNAstring(buffer, 1, d))
	    SET_STRING_ELT(ans, i, NA_STRING);
	else
	    SET_STRING_ELT(ans, i, insertString(buffer, d));
	break;
    case RAWSXP:
	if (isNAstring(buffer, 0, d))
	    RAW(ans)[i] = 0;
	else {
	    RAW(ans)[i] = strtoraw(buffer, &endp);
	    if (!isBlankString(endp))
		expected("a raw", buffer, d);
	}
	break;
    default:
	UNIMPLEMENTED_TYPE("extractItem", ans);
    }
}

SEXP R_plmd_model(SEXP Y, SEXP PsiCode, SEXP PsiK, SEXP Groups, SEXP Ngroups){


  SEXP R_return_value;
  SEXP R_weights;
  SEXP R_residuals;
  SEXP R_beta;
  SEXP R_SE;
  SEXP R_was_split;
  
  SEXP R_return_value_names;

  SEXP dim1;

  double *beta;
  double *residuals;
  double *weights;
  double *se;
  
  int *was_split;
  int *groups;

  double residSE;

  double *Ymat;

  double *X; /* Needed for SE */
  int X_cols, X_rows;


  int rows;
  int cols;

  int ngroups;
  
  int howmany_split =0;

  int i;
  
  PROTECT(dim1 = getAttrib(Y,R_DimSymbol));
  rows = INTEGER(dim1)[0];
  cols = INTEGER(dim1)[1];
  UNPROTECT(1);

  PROTECT(R_return_value = allocVector(VECSXP,5));
  
  /*
    Don't allocate R_beta/R_SE straight away, we won't know how much space
    these will actually need until finishing the PLM-d fitting procedure. 
    Instead we will just allocate those for which we currently know the size
  */
  
  
  PROTECT(R_weights = allocMatrix(REALSXP,rows,cols));
  PROTECT(R_residuals = allocMatrix(REALSXP,rows,cols));
  PROTECT(R_was_split = allocVector(INTSXP,rows));

  
  /* 0 - beta   (added below)
     1 - weights
     2 - residuals
     3 - standard errors  (added below)
     4 - R_was_split
  */
  SET_VECTOR_ELT(R_return_value,1,R_weights);
  SET_VECTOR_ELT(R_return_value,2,R_residuals);
  SET_VECTOR_ELT(R_return_value,4,R_was_split);



  UNPROTECT(3);

  residuals = NUMERIC_POINTER(R_residuals);
  weights = NUMERIC_POINTER(R_weights);
  was_split = INTEGER_POINTER(R_was_split);
  
  groups = INTEGER_POINTER(Groups);

  ngroups = INTEGER(Ngroups)[0];

  Ymat = NUMERIC_POINTER(Y);
  
  beta = Calloc(cols + rows*ngroups -1, double);   
  se = Calloc(cols + rows*ngroups -1, double);

  plmd_fit(Ymat, rows, cols, ngroups, groups, was_split, beta, residuals, weights, PsiFunc(asInteger(PsiCode)),asReal(PsiK), 20);


  for (i = 0; i < rows; i++){
    howmany_split+=was_split[i];
  }
  
  if (howmany_split > 0){
    PROTECT(R_beta = allocVector(REALSXP,rows + cols + howmany_split*(ngroups-1)));
    PROTECT(R_SE = allocVector(REALSXP,rows + cols + howmany_split*(ngroups-1)));
    
    X = plmd_get_design_matrix(rows, cols, ngroups, groups,was_split,&X_rows,&X_cols);


    rlm_compute_se(X,Ymat, X_rows, X_cols, beta, residuals, weights, se,(double *)NULL, &residSE, 2, PsiFunc(asInteger(PsiCode)),asReal(PsiK));
//.........这里部分代码省略.........

int layoutNRow(SEXP l) {
    return INTEGER(VECTOR_ELT(l, LAYOUT_NROW))[0];
}

//EXPORT int JpmcdsCdsoneUpfrontCharge(cdsone.c)
SEXP calcUpfrontTest
(SEXP baseDate_input,  /* (I) Value date  for zero curve       */
 SEXP types, /* "MMMMMSSSSSSSSS"*/
 SEXP rates, /* rates[14] = {1e-9, 1e-9, 1e-9, 1e-9, 1e-9, 1e-9, 1e-9,
		1e-9, 1e-9, 1e-9, 1e-9, 1e-9, 1e-9, 1e-9};/\* (I)
		Array of swap rates *\/ */
 SEXP expiries,
 SEXP mmDCC,          /* (I) DCC of MM instruments            */

 SEXP fixedSwapFreq,   /* (I) Fixed leg freqency/interval               */
 SEXP floatSwapFreq,   /* (I) Floating leg freqency/interval            */
 SEXP fixedSwapDCC,    /* (I) DCC of fixed leg                 */
 SEXP floatSwapDCC,    /* (I) DCC of floating leg              */
 SEXP badDayConvZC, //'M'  badDayConv for zero curve
 SEXP holidays,//'None'

 // input for upfront charge calculation
 SEXP todayDate_input, /*today: T (Where T = trade date)*/
 SEXP valueDate_input, /* value date: T+3 Business Days*/
 SEXP benchmarkDate_input,/* start date of benchmark CDS for internal
				     ** clean spread bootstrapping;
				     ** accrual Begin Date  */
 SEXP startDate_input,/* Accrual Begin Date */
 SEXP endDate_input,/*  Maturity (Fixed) */
 SEXP stepinDate_input,  /* T + 1*/
 
 SEXP dccCDS, 			/* accruedDcc */
 SEXP ivlCDS,
 SEXP stubCDS,
 SEXP badDayConvCDS,
 SEXP calendar,

 SEXP parSpread,
 SEXP couponRate,
 SEXP recoveryRate,
 SEXP isPriceClean_input,
 SEXP payAccruedOnDefault_input,
 SEXP notional) 

{
  //  static char routine[] = "JpmcdsCdsoneUpfrontCharge";

  // my vars
  int n;
  TDate baseDate, today, benchmarkDate, startDate, endDate, stepinDate,valueDate;
  int isPriceClean, payAccruedOnDefault;
  SEXP upfrontPayment;
  TCurve *discCurve = NULL;
  char* pt_types;
  char* pt_holidays;
  char* pt_mmDCC;
  char* pt_fixedSwapDCC;
  char* pt_floatSwapDCC;
  char* pt_fixedSwapFreq;
  char* pt_floatSwapFreq;
  char* pt_dccCDS;
  char* pt_ivlCDS;
  char* pt_stubCDS;
  char* pt_calendar;
  char* pt_badDayConvCDS;

  // new
  char *pt_badDayConvZC;
  double parSpread_for_upf, couponRate_for_upf, recoveryRate_for_upf, notional_for_upf;
  
  // function to consolidate R input to TDate
  baseDate_input = coerceVector(baseDate_input,INTSXP);
  baseDate = JpmcdsDate((long)INTEGER(baseDate_input)[0], 
			(long)INTEGER(baseDate_input)[1], 
			(long)INTEGER(baseDate_input)[2]);

  todayDate_input = coerceVector(todayDate_input,INTSXP);
  today = JpmcdsDate((long)INTEGER(todayDate_input)[0], 
		     (long)INTEGER(todayDate_input)[1], 
		     (long)INTEGER(todayDate_input)[2]);

  valueDate_input = coerceVector(valueDate_input,INTSXP);
  valueDate = JpmcdsDate((long)INTEGER(valueDate_input)[0], 
			 (long)INTEGER(valueDate_input)[1], 
			 (long)INTEGER(valueDate_input)[2]);

  benchmarkDate_input = coerceVector(benchmarkDate_input,INTSXP);
  benchmarkDate = JpmcdsDate((long)INTEGER(benchmarkDate_input)[0], 
			     (long)INTEGER(benchmarkDate_input)[1],
			     (long)INTEGER(benchmarkDate_input)[2]);

  startDate_input = coerceVector(startDate_input,INTSXP);
  startDate = JpmcdsDate((long)INTEGER(startDate_input)[0], 
			 (long)INTEGER(startDate_input)[1], 
			 (long)INTEGER(startDate_input)[2]);

  endDate_input = coerceVector(endDate_input,INTSXP);
  endDate = JpmcdsDate((long)INTEGER(endDate_input)[0],
		       (long)INTEGER(endDate_input)[1],
		       (long)INTEGER(endDate_input)[2]);

  stepinDate_input = coerceVector(stepinDate_input,INTSXP);
  stepinDate = JpmcdsDate((long)INTEGER(stepinDate_input)[0],
		       (long)INTEGER(stepinDate_input)[1],
//.........这里部分代码省略.........

SEXP glm_mcmcbas(SEXP Y, SEXP X, SEXP Roffset, SEXP Rweights, 
		 SEXP Rprobinit, SEXP Rmodeldim, 
		 SEXP modelprior, SEXP betaprior, SEXP Rbestmodel,SEXP plocal, 
		 SEXP BURNIN_Iterations,
		 SEXP family, SEXP Rcontrol,
		 SEXP Rupdate, SEXP Rlaplace) 
{
	int nProtected = 0;
	int nModels=LENGTH(Rmodeldim);

	SEXP ANS = PROTECT(allocVector(VECSXP, 17)); ++nProtected;
	SEXP ANS_names = PROTECT(allocVector(STRSXP, 17)); ++nProtected;
	SEXP Rprobs = PROTECT(duplicate(Rprobinit)); ++nProtected;
	SEXP MCMCprobs= PROTECT(duplicate(Rprobinit)); ++nProtected;
	SEXP R2 = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP shrinkage = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP modelspace = PROTECT(allocVector(VECSXP, nModels)); ++nProtected;
	SEXP modeldim =  PROTECT(duplicate(Rmodeldim)); ++nProtected;
	SEXP counts =  PROTECT(duplicate(Rmodeldim)); ++nProtected;
	SEXP beta = PROTECT(allocVector(VECSXP, nModels)); ++nProtected;
	SEXP se = PROTECT(allocVector(VECSXP, nModels)); ++nProtected;
	SEXP deviance = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP modelprobs = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP priorprobs = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP logmarg = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP sampleprobs = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP Q = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	SEXP Rintercept = PROTECT(allocVector(REALSXP, nModels)); ++nProtected;
	
	SEXP NumUnique = PROTECT(allocVector(INTSXP, 1)); ++nProtected;

	double *probs, MH=0.0, prior_m=1.0,shrinkage_m, logmargy, postold, postnew;
	int i, m, n, pmodel_old, *bestmodel;
	int mcurrent, n_sure;
	
	glmstptr *glmfamily;
	glmfamily = make_glmfamily_structure(family);

	betapriorptr *betapriorfamily;
	betapriorfamily = make_betaprior_structure(betaprior, family);


	//get dimsensions of all variables 
	int p = INTEGER(getAttrib(X,R_DimSymbol))[1];
	int k = LENGTH(modelprobs);
	int update = INTEGER(Rupdate)[0];
	double eps = DBL_EPSILON;

	struct Var *vars = (struct Var *) R_alloc(p, sizeof(struct Var)); // Info about the model variables. 
	probs =  REAL(Rprobs);
	n = sortvars(vars, probs, p); 
	for (i =n; i <p; i++) REAL(MCMCprobs)[vars[i].index] = probs[vars[i].index];
	for (i =0; i <n; i++) REAL(MCMCprobs)[vars[i].index] = 0.0;

	// fill in the sure things 
	int *model = ivecalloc(p);
	for (i = n, n_sure = 0; i < p; i++)  {
		model[vars[i].index] = (int) vars[i].prob;
		if (model[vars[i].index] == 1) ++n_sure;
	}

	GetRNGstate();

	NODEPTR tree, branch;
	tree = make_node(-1.0);
	//  Rprintf("For m=0, Initialize Tree with initial Model\n");  

	m = 0;
	bestmodel = INTEGER(Rbestmodel);
	INTEGER(modeldim)[m] = n_sure;

	// Rprintf("Create Tree\n"); 
	branch = tree;
	CreateTree(branch, vars, bestmodel, model, n, m, modeldim);
	int pmodel = INTEGER(modeldim)[m];
	SEXP Rmodel_m =	PROTECT(allocVector(INTSXP,pmodel));
	GetModel_m(Rmodel_m, model, p);
	//evaluate logmargy and shrinkage
	SEXP glm_fit = PROTECT(glm_FitModel(X, Y, Rmodel_m, Roffset, Rweights,
					    glmfamily, Rcontrol, Rlaplace,
					    betapriorfamily));	
	prior_m  = compute_prior_probs(model,pmodel,p, modelprior);

	logmargy = REAL(getListElement(getListElement(glm_fit, "lpy"),"lpY"))[0];
	shrinkage_m = REAL(getListElement(getListElement(glm_fit, "lpy"),	
						  "shrinkage"))[0];
	SetModel2(logmargy, shrinkage_m, prior_m, sampleprobs, logmarg, shrinkage, priorprobs, m);
	SetModel1(glm_fit, Rmodel_m, beta, se, modelspace, deviance, R2,
		  Q, Rintercept, m);
	UNPROTECT(2);

	int nUnique=0, newmodel=0;
	double *real_model = vecalloc(n);
	int *modelold = ivecalloc(p);
	int old_loc = 0;
	int new_loc;
	pmodel_old = pmodel;
	nUnique=1;
	INTEGER(counts)[0] = 0;
	postold =  REAL(logmarg)[m] + log(REAL(priorprobs)[m]);
//.........这里部分代码省略.........

SEXP BATCHgdwm_loop_MLPnet (SEXP origNet, SEXP Ptrans, SEXP Ttrans, SEXP nepochs, SEXP rho, SEXP thread_number ) {
   //The only difference between wm and without it is the weight update (and the place the values are stored, one is batchgd the other is batchgdwm)
   SEXP net;
   SEXP R_fcall, args, arg1, arg2, arg3;

   PROTECT(net=duplicate(origNet));
   int* Ptransdim = INTEGER(coerceVector(getAttrib(Ptrans, R_DimSymbol), INTSXP));
   int* Ttransdim = INTEGER(coerceVector(getAttrib(Ttrans, R_DimSymbol), INTSXP));
   int n_epochs  = INTEGER(nepochs)[0];
   struct AMOREnet* ptnet = copynet_RC(net);
   struct AMOREneuron** neurons = ptnet->neurons;

   /////////////////////////////////////////////////////////////////////////
   //Convert input and target to double only once (and instead of copying it every time, just change the pointers)
   //Different rows for easy switching pointers
   double*  input_data  = REAL(Ptrans);
   double*  target_data = REAL(Ttrans);
   double** inputs  = (double**) R_alloc(Ptransdim[1],sizeof(double*)); //This is an 'Index'
   double** targets = (double**) R_alloc(Ptransdim[1],sizeof(double*)); //This is an 'Index'

   for (int fila=0; fila < Ptransdim[1]; fila++) {
      inputs[fila]  = &input_data [fila*Ptransdim[0]];
      targets[fila] = &target_data[fila*Ttransdim[0]];
   }
   /////////////////////////////////////////////////////////////////////////

   /////////////////////////////////////////////////////////////////////////
   // Thread number calculation
   int n_threads = 1;
#ifdef _OPENMP
   {
      int max_threads = omp_get_max_threads();
      int given_threads = 0;

      if (isInteger(thread_number))
        given_threads = INTEGER(thread_number)[0];
      else if (isNumeric(thread_number))
        given_threads = floor(REAL(thread_number)[0]);

      if (given_threads <1) //I HAVE THE POWER TO SCHEDULE!
        if(max_threads  >1)
          n_threads = max_threads-1; //Leave a CPU free
        else
          n_threads = 1;
      else if (given_threads > max_threads)
        n_threads = max_threads;
      else
        n_threads = given_threads;

      if (neurons[0]->actf == CUSTOM_ACTF) //OMP + R custom functions = bad idea
        n_threads = 1;
      else if ((ptnet->deltaE.name != LMLS_NAME) && (ptnet->deltaE.name != LMS_NAME))
        n_threads = 1;

      //printf("Using %i threads from a max of %i.\n",n_threads ,max_threads);
   }
#endif
   /////////////////////////////////////////////////////////////////////////

   /////////////////////////////////////////////////////////////////////////
   //Contribution (who is to blame) : Parallelization done by Jose Maria
   //Memory allocation for running different threads in parallel:
   // Each thread will have his own pool of memory to handle the two kinds of temp vars:
   //   Vars used only inside the forwards/backwards (v0, v1 and method_delta)
   //     These vars will be initialized and read only by each thread
   //   Vars that hold the information on how much the weights and the bias should change 
   //     These vars will be initialized by each thread, then accumulated and read by the master thread when the batch is finished
   int n_neurons = ptnet->last_neuron+1;
   //Temp values, internal in each iteration
   double **  v0s                 = (double** ) R_alloc(n_threads,sizeof(double* )); //This is an 'Index'
   double **  v1s                 = (double** ) R_alloc(n_threads,sizeof(double* )); //This is an 'Index'
   double **  method_deltas       = (double** ) R_alloc(n_threads,sizeof(double* )); //This is an 'Index'
   //Accumulated values
   double **  method_deltas_bias  = (double** ) R_alloc(n_threads,sizeof(double* )); //This is an 'Index'
   double *** method_sums_delta_x = (double***) R_alloc(n_threads,sizeof(double**)); //This is an 'Index'

   for(int id_thread=0; id_thread<n_threads;id_thread++){
      double* chunk = (double*) R_alloc(4*n_neurons,sizeof(double)); //Actual chunk of memory for each thread, trying to avoid R_alloc calls
      //Advantages: Good proximity reference in cache for values of the same thread, and since it has at least 2 neurons
      // (Who would have a NNetwork with less than 2 neurons?), chunks are larger than 64 bytes (i7 L2 cache block size?)
      v0s               [id_thread] =  chunk             ;  
      v1s               [id_thread] = &chunk[  n_neurons];
      method_deltas     [id_thread] = &chunk[2*n_neurons];
      method_deltas_bias[id_thread] = &chunk[3*n_neurons];
      
      method_sums_delta_x[id_thread] = (double**) R_alloc(n_neurons,sizeof(double*)); //This is an 'Index'
      for(int i=0; i<n_neurons; i++) //Different weigth number for each layer, TODO: R_alloc each layer instead of each neuron
         method_sums_delta_x[id_thread][i] = (double*) R_alloc(neurons[i]->last_input_link+1,sizeof(double));
   }
   /////////////////////////////////////////////////////////////////////////

   /////////////////////////////////////////////////////////////////////////
   //Consistency (leave pnet as if the function had worked with their values instead of external ones)
   // R_alloc should handle freeing the memory, so it's not needed to free the previously allocated memory to avoid memory leaks
   // Changing pointer instead of copying data
   ptnet->input  = inputs[Ptransdim[1]-1];
   ptnet->target = targets[Ptransdim[1]-1];
   /////////////////////////////////////////////////////////////////////////
   
   /////////////////////////////////////////////////////////////////////////
//.........这里部分代码省略.........

/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 
 * Purpose: 
 *          
 *          
 *          
 *
 *          
*/
SEXP write_volume(SEXP filename, SEXP nDimensions, 
								SEXP dimLengths, 
								SEXP dimStarts, 
								SEXP dimSteps, 
								SEXP volumeDataType, 
								SEXP volumeRange, 
								SEXP hSlab) {

	mihandle_t			minc_volume;
	midimhandle_t		dim[MI2_MAX_VAR_DIMS];
	mivolumeprops_t		volume_properties;
	
	int				result;
	int				ndx;
	int				no_dimensions;
	int				dim_lengths[MI2_MAX_VAR_DIMS];
	double			dim_starts[MI2_MAX_VAR_DIMS];
	double			dim_steps[MI2_MAX_VAR_DIMS];
	double			volume_range_min, volume_range_max;
	int				volume_data_type;
	
	// pointer to the volume data
	double			*hSlab_p;
	misize_t	hSlab_start[MI2_MAX_VAR_DIMS];
	misize_t	hSlab_count[MI2_MAX_VAR_DIMS];


	// start ...
	if ( R_DEBUG_mincIO ) Rprintf("write_volume: start ...\n");

	/* init number of dimensions and their respective sizes
	  ... yes, I could get this myself, but it's best set in the calling R code */
	no_dimensions = INTEGER(nDimensions)[0];
	volume_data_type = INTEGER(volumeDataType)[0];
	volume_range_min = REAL(volumeRange)[0];
	volume_range_max = REAL(volumeRange)[1];
	
	// init volume data pointer
	hSlab_p = REAL(hSlab);
	
	
	// init lenghts, steps, and starts
	for (ndx=0; ndx < no_dimensions; ++ndx) {
		dim_lengths[ndx]  = INTEGER(dimLengths)[ndx];
		hSlab_count[ndx]  = INTEGER(dimLengths)[ndx];
		dim_starts[ndx]  = REAL(dimStarts)[ndx];
		dim_steps[ndx]  = REAL(dimSteps)[ndx];
	}


	// set the properties for the new output volume
	// ... no compression, no chunking, no multi-resolution, ... nothin fancy
	result = minew_volume_props(&volume_properties);
	if (result != MI_NOERROR) {
		error("write_volume:minew_volume_props: Error setting output volume properties: %s.\n", CHAR(STRING_ELT(filename, 0)));
	}
	result = miset_props_compression_type(volume_properties, MI_COMPRESS_NONE);
	if (result != MI_NOERROR) {
		error("write_volume:miset_props_compression_type: Error setting output volume properties: %s.\n", CHAR(STRING_ELT(filename, 0)));
	}
	result = miset_props_multi_resolution(volume_properties, FALSE, 1);
	if (result != MI_NOERROR) {
		error("write_volume:miset_props_multi_resolution: Error setting output volume properties: %s.\n", CHAR(STRING_ELT(filename, 0)));
	}

	/*	create the 3 output dimensions in the order Z, Y, X, as the volume
		is stored in this order */
	// z-dim
	result = micreate_dimension("zspace", 
								MI_DIMCLASS_SPATIAL,
								MI_DIMATTR_REGULARLY_SAMPLED, 
								dim_lengths[0],
								&dim[0]);
	//
	if (result != MI_NOERROR) {
		error("write_volume: Error initializing the dimension struct for %s dimension.\n", "zspace");
	}
	
	// y-dim
	result = micreate_dimension("yspace", 
								MI_DIMCLASS_SPATIAL,
								MI_DIMATTR_REGULARLY_SAMPLED, 
								dim_lengths[1],
								&dim[1]);
	//
	if (result != MI_NOERROR) {
		error("write_volume: Error initializing the dimension struct for %s dimension.\n", "yspace");
	}
	
	// x-dim
	result = micreate_dimension("xspace", 
								MI_DIMCLASS_SPATIAL,
//.........这里部分代码省略.........

QueryData genMDDevices(QueryContext& context) {
  QueryData results;
  MDStat mds;
  MD md;
  std::vector<std::string> lines;

  getLines(lines);

  md.parseMDStat(lines, mds);
  for (const auto& device : mds.devices) {
    std::string path(md.getPathByDevName(device.name));
    if (path.empty()) {
      LOG(ERROR) << "Could not get file path for " << device.name;
      return results;
    }

    mdu_array_info_t array;
    if (!md.getArrayInfo(path, array)) {
      return results;
    }

    Row r;
    r["device_name"] = device.name;
    r["status"] = device.status;
    r["raid_level"] = INTEGER(array.level);
    r["size"] = BIGINT(device.usableSize);
    r["chunk_size"] = BIGINT(array.chunk_size);
    r["raid_disks"] = INTEGER(array.raid_disks);
    r["nr_raid_disks"] = INTEGER(array.nr_disks);
    r["working_disks"] = INTEGER(array.working_disks);
    r["active_disks"] = INTEGER(array.active_disks);
    r["failed_disks"] = INTEGER(array.failed_disks);
    r["spare_disks"] = INTEGER(array.spare_disks);

    r["superblock_state"] = getSuperBlkStateStr(array.state);
    r["superblock_version"] = md.getSuperblkVersion(device.name);
    r["superblock_update_time"] = BIGINT(array.utime);

    if (!device.recovery.progress.empty()) {
      r["recovery_progress"] = device.recovery.progress;
      r["recovery_finish"] = device.recovery.finish;
      r["recovery_speed"] = device.recovery.speed;
    }

    if (!device.resync.progress.empty()) {
      r["resync_progress"] = device.resync.progress;
      r["resync_finish"] = device.resync.finish;
      r["resync_speed"] = device.resync.speed;
    }

    if (!device.reshape.progress.empty()) {
      r["reshape_progress"] = device.reshape.progress;
      r["reshape_finish"] = device.reshape.finish;
      r["reshape_speed"] = device.reshape.speed;
    }

    if (!device.checkArray.progress.empty()) {
      r["check_array_progress"] = device.checkArray.progress;
      r["check_array_finish"] = device.checkArray.finish;
      r["check_array_speed"] = device.checkArray.speed;
    }

    if (!device.bitmap.onMem.empty()) {
      r["bitmap_on_mem"] = device.bitmap.onMem;
      r["bitmap_chunk_size"] = device.bitmap.chunkSize;
      r["bitmap_external_file"] = device.bitmap.externalFile;
    }

    r["other"] = device.other;
    r["unused_devices"] = mds.unused;

    results.push_back(r);
  }

  return results;
}

//extern "C"
SEXP mc_irf_var(SEXP varobj, SEXP nsteps, SEXP draws)
{
  int m, p, dr=INTEGER(draws)[0], ns=INTEGER(nsteps)[0], T, df, i;
  SEXP AR, Y, Bhat, XR, prior, hstar, meanS, output;

  // Get # vars/lags/steps/draws/T/df
  PROTECT(AR = listElt(varobj, "ar.coefs"));
  PROTECT(Y = listElt(varobj, "Y"));
  m = INTEGER(getAttrib(AR, R_DimSymbol))[0]; //#vars
  p = INTEGER(getAttrib(AR, R_DimSymbol))[2]; //#lags
  T = nrows(Y); df = T - m*p - m - 1;
  UNPROTECT(2);

  // Put coefficients from varobj$Bhat in Bcoefs vector (m^2*p, 1)
  PROTECT(Bhat = coerceVector(listElt(varobj, "Bhat"), REALSXP));
  Matrix bcoefs = R2Cmat(Bhat, m*p, m);
  bcoefs = bcoefs.AsColumn();
  UNPROTECT(1);

  // Define X(T x m*p) subset of varobj$X and XXinv as solve(X'X)
  PROTECT(XR = coerceVector(listElt(varobj,"X"),REALSXP));
  Matrix X = R2Cmat(XR, T, m*p), XXinv;
  UNPROTECT(1);

  // Get the correct moment matrix
  PROTECT(prior = listElt(varobj,"prior"));
  if(!isNull(prior)){
    PROTECT(hstar = coerceVector(listElt(varobj,"hstar"),REALSXP));
    XXinv = R2Cmat(hstar, m*p, m*p).i();
    UNPROTECT(1);
  }
  else { XXinv = (X.t()*X).i(); }
  UNPROTECT(1);

  // Get the transpose of the Cholesky decomp of XXinv
  SymmetricMatrix XXinvSym; XXinvSym << XXinv;
  XXinv = Cholesky(XXinvSym);

  // Cholesky of covariance
  PROTECT(meanS = coerceVector(listElt(varobj,"mean.S"),REALSXP));
  SymmetricMatrix meanSSym; meanSSym << R2Cmat(meanS, m, m);
  Matrix Sigmat = Cholesky(meanSSym);
  UNPROTECT(1);

  // Matricies needed for the loop
  ColumnVector bvec; bvec=0.0;
  Matrix sqrtwish, impulse(dr,m*m*ns); impulse = 0.0;
  SymmetricMatrix sigmadraw; sigmadraw = 0.0;
  IdentityMatrix I(m);

  GetRNGstate();
  // Main Loop
  for (i=1; i<=dr; i++){
    // Wishart/Beta draws
    sigmadraw << Sigmat*(T*rwish(I,df).i())*Sigmat.t();
    sqrtwish = Cholesky(sigmadraw);
    bvec = bcoefs+KP(sqrtwish, XXinv)*rnorms(m*m*p);

    // IRF computation
    impulse.Row(i) = irf_var_from_beta(sqrtwish, bvec, ns).t();
    if (!(i%1000)){ Rprintf("Monte Carlo IRF Iteration = %d\n",i); }
  } // end main loop
  PutRNGstate();

  int dims[]={dr,ns,m*m};
  PROTECT(output = C2R3D(impulse,dims));
  setclass(output,"mc.irf.VAR");
  UNPROTECT(1);
  return output;
}