   float *dev_proj)
 {
   // Copy the size of inputs into constant memory
-  cudaMemcpyToSymbol(c_projSize, proj_size, sizeof(int3));
+  cudaMemcpyToSymbol(c_proj_size, proj_size, sizeof(int3));
   cudaMemcpyToSymbol(c_vol_size, vol_size, sizeof(int3));
 
   // Copy the projection matrices into constant memory

  float3 ip;
  float  voxel_data = 0;

  for (unsigned int proj = 0; proj<c_projSize.z; proj++)
    {
    // matrix multiply
    ip = matrix_multiply(make_float3(i,j,k), &(c_matrices[12*proj]));


  // Place it into the volume
  dev_vol_out[vol_idx] = dev_vol_in[vol_idx] + voxel_data;
}


__global__
void kernel_optim(float *dev_vol_in, float *dev_vol_out)
{
  // Allocate a few shared buffers
  extern __shared__ float shared[];
  float* accumulators = &shared[0];

  // Compute the starting position of the shared memory available to the current thread
  unsigned int startIndex = (threadIdx.x + (threadIdx.y * ( blockDim.x + threadIdx.z * blockDim.y) ) ) * 8;

  unsigned int i = blockIdx.x * blockDim.x + threadIdx.x;
  unsigned int j = (blockIdx.y * blockDim.y + threadIdx.y)*8;
  unsigned int k = blockIdx.z * blockDim.z + threadIdx.z;

  if (i >= c_vol_size.x || j >= c_vol_size.y || k >= c_vol_size.z)
    {
    return;
    }
  unsigned int range_j = max(c_vol_size.y - j, 8);

  // Index row major into the volume
  long int vol_idx = i + (j + k*c_vol_size.y)*(c_vol_size.x);

  float3 ip;

  // Initialize the shared memory with zeros
  for(unsigned int l=0; l<range_j; l++)
    accumulators[startIndex + l]=0;

  // Process all the projections in the stack
  for (unsigned int proj = 0; proj<c_proj_size.z; proj++)
    {
    // matrix multiply
    ip = matrix_multiply(make_float3(i,j,k), &(c_matrices[12*proj]));

    // Change coordinate systems
    ip.z = 1 / ip.z;
    ip.x = ip.x * ip.z;
    ip.y = ip.y * ip.z;
    float dx = c_matrices[12*proj +  1]*ip.z;
    float dy = c_matrices[12*proj +  5]*ip.z;

    // Walk voxels and accumulate their back projection in shared memory
    for(unsigned int l=0; l<range_j; l++)
      {
      accumulators[startIndex + l] += tex2DLayered(tex_proj, ip.x, ip.y, proj);
      ip.x+=dx;
      ip.y+=dy;
      }
    }

  // Write the accumulated results into global memory
  for(unsigned int l=0; l<range_j; l++)
    {
    dev_vol_out[vol_idx] = dev_vol_in[vol_idx] + accumulators[startIndex + l];
    vol_idx+=c_vol_size.x;
    }
}

//_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
