     t30 = totals[threadIdx.x][30];
     t31 = totals[threadIdx.x][31];
 
+    // The first feature index is actually the floating point threshold
     fthresh = *((float *)&t0);
 
     // Clear totals for each tree
     totals[threadIdx.x][0] = 0;
   
     // Now read the column and update totals
-    for (tid = 0; tid + threadIdx.x < nrows; tid += blockDim.x) {
-      vv = feats[tid + threadIdx.x + bd * nrows];
-      imin = tid;
-      imax = tid + blockDim.x;
-      tx = min(31, max(0, t1 - imin)); vtmp = __shfl(vv, tx); if (t1 >= imin && t1 < imax) totals[threadIdx.x][0] += vtmp;
+    for (tid = 0; tid < nrows; tid += blockDim.x) {
+      if (tid + threadIdx.x < nrows) {
+        vv = feats[tid + threadIdx.x + bd * nrows];
+      }
+      //      imin = tid;
+      //      imax = tid + blockDim.x;
+      // Check each feature index in turn, and see if it matches an input feature. Skip t0 which is actually the thresholds.
+      /*      tx = min(31, max(0, t1 - imin)); vtmp = __shfl(vv, tx); if (t1 >= imin && t1 < imax) totals[threadIdx.x][0] += vtmp;
       tx = min(31, max(0, t2 - imin)); vtmp = __shfl(vv, tx); if (t2 >= imin && t2 < imax) totals[threadIdx.x][0] += vtmp;
       tx = min(31, max(0, t3 - imin)); vtmp = __shfl(vv, tx); if (t3 >= imin && t3 < imax) totals[threadIdx.x][0] += vtmp;
       tx = min(31, max(0, t4 - imin)); vtmp = __shfl(vv, tx); if (t4 >= imin && t4 < imax) totals[threadIdx.x][0] += vtmp;

      tx = min(31, max(0, t28 - imin)); vtmp = __shfl(vv, tx); if (t28 >= imin && t28 < imax) totals[threadIdx.x][0] += vtmp;
      tx = min(31, max(0, t29 - imin)); vtmp = __shfl(vv, tx); if (t29 >= imin && t29 < imax) totals[threadIdx.x][0] += vtmp;
      tx = min(31, max(0, t30 - imin)); vtmp = __shfl(vv, tx); if (t30 >= imin && t30 < imax) totals[threadIdx.x][0] += vtmp;
      tx = min(31, max(0, t31 - imin)); vtmp = __shfl(vv, tx); if (t31 >= imin && t31 < imax) totals[threadIdx.x][0] += vtmp;
    }

    newt = 2 * tp;
    if (totals[threadIdx.x][0] > fthresh) {
      newt++;

    t13 = totals[threadIdx.x][13];
    t14 = totals[threadIdx.x][14];
    t15 = totals[threadIdx.x][15];

    fthresh = *((float *)&t0);

