               float *input_p = input + ((long)table[ii]-1)*input_h*input_w 
                 + yy*stride_h*input_w + xx*stride_w;
               float *output_p = output + oo*output_h*output_w + yy*output_w + xx;
-              float *kernel_p = kernel + ((long)table[ii + 1]-1) *kernel_w*kernel_h + koffset;
+              float *kernel_p = shared_kernel 
+                + ((long)table[ii + 1]-1) *kernel_w*kernel_h + koffset;
+              float sum = 0;
+              if (swapkernel) {
+#pragma unroll
+                for(ky = 0; ky < T_kernel_h; ky++) {
+#pragma unroll
+                  for(kx = 0; kx < T_kernel_w; kx++) {
+                    sum += input_p[kx]*(*kernel_p--);
+                  }
+                  input_p += input_w;
+                }
+              } else {
+#pragma unroll
+                for(ky = 0; ky < T_kernel_h; ky++) {
+#pragma unroll
+                  for(kx = 0; kx < T_kernel_w; kx++) {
+                    sum += input_p[kx]*(*kernel_p++);
+                  }
+                  input_p += input_w;
+                }
+              }
+              *output_p += sum;
+            }
+          }
+        }
+      }
+    } else {
+      // default convolution loop
+      for(oo = oo_start; oo < oo_end; oo++) {
+        for (ii = table_start; ii < table_end; ii++) {
+          for(yy = yy_start; yy < yy_end; yy+=yy_step) {
+            for(xx = xx_start; xx < xx_end; xx+=xx_step) {
+              // Dot product in two dims (between input image and the mask)
+              float *input_p = input + ((long)table[ii]-1)*input_h*input_w 
+                + yy*stride_h*input_w + xx*stride_w;
+              float *output_p = output + oo*output_h*output_w + yy*output_w 
+                + xx;
+              float *kernel_p = shared_kernel 
+                + (((long)table[ii]-1) % fanin) * kernel_w * kernel_h + koffset;
               float sum = 0;
               if (swapkernel) {
                 for(ky = 0; ky < kernel_h; ky++) {


///////////////////////////////////
///// ConvolutionMap
  /*
   * Description:
   *   base conv2D routine: 3D input, 3D output, 4D kernel
   *
   *   - all chunks of data should be contiguous
   *   - the swapkernel flag can be used to generate a conv2 instead of xcorr2
   *   - the templated kernel size is useful to generate code that's 2x faster
   *     but can be set to 0 to allow arbitrary kernel sizes
   *   ---- the table should have the first dim with the outputs, each output 
   *   ---- should have a fanin set of inputs contiguously
   */
  template <bool swapkernel, int T_kernel_h, int T_kernel_w>
  __global__ void conv2mapgeneric(float *input, float *kernel, float *output,
                                  int input_n, int input_h, int input_w,
                                  int kernel_n, int kernel_h, int kernel_w,
                                  int stride_w, int stride_h,
                                  float *table, int fanin)
  {
    // output dimensions
    int output_h = (input_h - kernel_h) / stride_h + 1;
    int output_w = (input_w - kernel_w) / stride_w + 1;

    // xcorr or conv
    int koffset = swapkernel ? kernel_w*kernel_h-1 : 0;

    // nb outputs
    // int output_n = kernel_n / fanin;

    // generate offsets according to block/thread ids
    int xx_start = threadIdx.x;
    int xx_end = output_w;
    int xx_step = blockDim.x;

    int yy_start = blockDim.y*blockIdx.y + threadIdx.y;
    int yy_end = output_h;
    int yy_step = blockDim.y*gridDim.y;

    int oo_start = blockIdx.x;
    int oo_end = oo_start+1;

    int table_start = blockIdx.x * (fanin * 2);
    int table_end = table_start + (fanin * 2);

    // nb threads, unique thread id
    int tid = blockDim.x*blockDim.y*threadIdx.z 
      + blockDim.x * threadIdx.y + threadIdx.x;
    int nthreads = blockDim.x * blockDim.y * blockDim.z;

    // iterators
    int oo, ii, xx, yy, kx, ky, kk;

    // do the kernels fit in shared mem ?
    if (kernel_w*kernel_h*kernel_n <= CUDA_SHARED_MEM_SIZE) { 
      // put the kernel in shared memory
      __shared__ float shared_kernel[CUDA_SHARED_MEM_SIZE];

      // first thread of each block does the copy
      for (kk = tid; kk < kernel_w*kernel_h*kernel_n; kk += nthreads) {
        shared_kernel[kk] = kernel[kk];
      }
      __syncthreads();

      // templated kernel size
      if ((T_kernel_w > 0) && (T_kernel_h > 0)) {
        // unrolled convolution loop
        for(oo = oo_start; oo < oo_end; oo++) {
          for (ii = table_start; ii < table_end; ii = ii + 2) {
            for(yy = yy_start; yy < yy_end; yy+=yy_step) {
              for(xx = xx_start; xx < xx_end; xx+=xx_step) {
                // Dot product in two dimensions... (between input image and the mask)
                float *input_p = input + ((long)table[ii]-1)*input_h*input_w 
                  + yy*stride_h*input_w + xx*stride_w;
                float *output_p = output + oo*output_h*output_w + yy*output_w + xx;
                float *kernel_p = shared_kernel 
                  + ((long)table[ii + 1]-1) *kernel_w*kernel_h + koffset;
                float sum = 0;
                if (swapkernel) {
#pragma unroll
                  for(ky = 0; ky < T_kernel_h; ky++) {
#pragma unroll
                    for(kx = 0; kx < T_kernel_w; kx++) {
                      sum += input_p[kx]*(*kernel_p--);
                    }
                    input_p += input_w;
                  }
                } else {
#pragma unroll
                  for(ky = 0; ky < T_kernel_h; ky++) {
#pragma unroll
                    for(kx = 0; kx < T_kernel_w; kx++) {
                      sum += input_p[kx]*(*kernel_p++);
                    }
                    input_p += input_w;
                  }
                }
                *output_p += sum;
              }
            }
          }
        }
      } else {
        // default convolution loop
        for(oo = oo_start; oo < oo_end; oo++) {
          for (ii = table_start; ii < table_end; ii++) {
            for(yy = yy_start; yy < yy_end; yy+=yy_step) {
              for(xx = xx_start; xx < xx_end; xx+=xx_step) {
                // Dot product in two dims (between input image and the mask)
                float *input_p = input + ((long)table[ii]-1)*input_h*input_w 
                  + yy*stride_h*input_w + xx*stride_w;
                float *output_p = output + oo*output_h*output_w + yy*output_w 
                  + xx;
                float *kernel_p = shared_kernel 
                  + (((long)table[ii]-1) % fanin) * kernel_w * kernel_h + koffset;
                float sum = 0;
                if (swapkernel) {
                  for(ky = 0; ky < kernel_h; ky++) {
#pragma unroll 5
                    for(kx = 0; kx < kernel_w; kx++) {
                      sum += input_p[kx]*(*kernel_p--);
                    }
                    input_p += input_w;
                  }
                } else {
                  for(ky = 0; ky < kernel_h; ky++) {
#pragma unroll 5
                    for(kx = 0; kx < kernel_w; kx++) {
                      sum += input_p[kx]*(*kernel_p++);
                    }
                    input_p += input_w;
                  }
                }
                *output_p += sum;
              }
            }
          }
        }
      }

    } else { // not enough shared mem for kernels, simply stream them

      // convolution loop
      for(oo = oo_start; oo < oo_end; oo++) {
        for (ii = table_start; ii < table_end; ii = ii + 2) {
          for(yy = yy_start; yy < yy_end; yy+=yy_step) {

                                                           stride_y, 
                                                           table_data, 
                                                           fanin);
  else
    conv2mapgeneric <false, 0 , 0> <<<blocks, threads>>> (input_data, 
                                                          kernel_data, 
                                                          output_data,
                                                          nInputPlane, 
                                                          nInputRows, 
                                                          nInputCols,
                                                          nOutputPlane*fanin, 
                                                          nKernelRows, 
                                                          nKernelCols,
                                                          stride_x, 
                                                          stride_y, 
                                                          table_data, 
                                                          fanin);
  // sync & clean
  cudaDeviceSynchronize();

  THCudaTensor_free(input);
  THCudaTensor_free(kernel);
  THCudaTensor_free(table);

  // check for errors
  cudaError_t err = cudaGetLastError();
