title: cuda convolution mapping
body: <p>i'm trying to copy for each block of threads a patch of image and relative apron to shared memory.</p>
  
  <p>After my data are copyied(i used a <strong>matrix</strong>) to shared memory, i want a relations that map the center of the mask in shared memory that i consider for convolution and the center of the mask in the image buffer. </p>
  
  <p>I want that because if i try to do convolution of image seems that the center of the mask in shared memory doesn't correspond to the center in the image buffer stored in <strong>global memory</strong>.</p>
  
  <p>In the code below i write an example of simple image <strong>black and white erosion algorithm</strong> , when i put the result of a <strong>convolution</strong> to the output image seems that the center not corresponds. </p>
  
  <p>i write my sample using a <strong>512x512</strong>px image</p>
  
  <p>my settings are:</p>
  
  <pre><code>//block and grid size
  dim3 block(16,16);
  dim3 grid(512/(block.x),512/(block.y),1);
  </code></pre>
  
  <p>this is my kernel:</p>
  
  <pre><code>#define STREL_SIZE 5
  
  #define TILE_W 16
  #define TILE_H 16
  
  #define R (STREL_SIZE/2)
  
  //size of tile image + apron
  #define BLOCK_W (TILE_W+(2*R))
  #define BLOCK_H (TILE_H+(2*R))
  
  
   __global__ void erode_multiple_img_SM_v2(unsigned char * buffer_in,
                              unsigned char * buffer_out,
                              int w,int h ){
  
      // Data cache: threadIdx.x , threadIdx.y
      __shared__ unsigned char data[TILE_W +STREL_SIZE  ][TILE_H +STREL_SIZE ];
  
  
       int col = blockIdx.x * blockDim.x + threadIdx.x;
       int row = blockIdx.y * blockDim.y + threadIdx.y;
  
       // global mem address of this thread
       int gLoc =  row*w +col;
  
  
       int x, y;  // image based coordinate
  
  
  
       if((col&lt;w)&amp;&amp;(row&lt;h)) {
           data[threadIdx.x][threadIdx.y]=buffer_in[gLoc];
  
       if (threadIdx.y &gt; (h-STREL_SIZE))
            data[threadIdx.x][threadIdx.y + STREL_SIZE]=buffer_in[gLoc + STREL_SIZE];
  
       if (threadIdx.x &gt;(w-STREL_SIZE))
            data[threadIdx.x + STREL_SIZE][threadIdx.y]=buffer_in[gLoc+STREL_SIZE];
  
       if ((threadIdx.x &gt;(w-STREL_SIZE)) &amp;&amp; (threadIdx.y &gt; (h-STREL_SIZE)))
            data[threadIdx.x+STREL_SIZE][threadIdx.y+STREL_SIZE] =     buffer_in[gLoc+2*STREL_SIZE];
  
       //wait for all threads to finish read
       __syncthreads();
  
        unsigned char min_value = 255;
        for(x=0;x&lt;STREL_SIZE;x++){
            for(y=0;y&lt;STREL_SIZE;y++){
                min_value = min( (data[threadIdx.x+x][threadIdx.y+y]) , min_value);
                }
  
            }
        buffer_out[gLoc]= min_value;
        }
  }
  </code></pre>
  
  <p>my input image:</p>
  
  <p><a href="https://i.stack.imgur.com/r7619.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/r7619.jpg" alt="input"></a></p>
  
  <p>my the output of the kernel is:</p>
  
  <p><a href="https://i.stack.imgur.com/STw86.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/STw86.jpg" alt="output"></a></p>
  
  <p>where <strong>w</strong> is the width of image,and is equal <strong>512</strong>,
  where <strong>h</strong> is the height of image,and is equal <strong>512</strong>.</p>
  
  <p>i call the kernel with:</p>
  
  <pre><code> erode_multiple_img_SM&lt;&lt;&lt;grid,block&gt;&gt;&gt;(dimage_src,dimage_dst,512,512);
  </code></pre>
  
  <p>the <strong>dimage_src</strong> is the input image an array buffer not a matrix, and <strong>dimage_dst</strong> is the output image a buffer.</p>
  
  <p>each buffer have the size of <strong>nElem * nImg * sizeof(unsigned char)</strong> where <strong>nElem=512*512</strong> is the size of the buffer and <strong>nImg</strong> is the number of image that i want processing in my case is equal to <strong>1</strong>. 
  where i'm wrong?</p>
  
  <p><strong>CODE UPDATE:</strong></p>
  
  <pre><code>__global__ void erode_multiple_img_SM_v2(unsigned char * buffer_in,
                              unsigned char * buffer_out,
                              int w,int h ){
  
  // Data cache: threadIdx.x , threadIdx.y
  __shared__ unsigned char data[TILE_W + STREL_SIZE-1 ][TILE_H + STREL_SIZE-1 ];
  
  // global mem address of this thread
  int col = blockIdx.x * blockDim.x + threadIdx.x;
  int row = blockIdx.y * blockDim.y + threadIdx.y;
  
  
  
  int gLoc =  row*w +col;
  
  
  
  // each threads loads four values from global memory into shared mem
  int x, y;   // image based coordinate
  
  
  
  if((col&lt;w)&amp;&amp;(row&lt;h)) {
  
      data[threadIdx.x][threadIdx.y] = buffer_in[gLoc];
  
       if (threadIdx.y &gt; (TILE_H-STREL_SIZE+1))
            data[threadIdx.x][threadIdx.y + STREL_SIZE-1]=buffer_in[(row + STREL_SIZE-1)*w + col];
  
       if (threadIdx.x &gt; (TILE_W-STREL_SIZE+1))
             data[threadIdx.x + STREL_SIZE-1][threadIdx.y] = buffer_in[row*w+col + STREL_SIZE-1];
  
       if ((threadIdx.x &gt; (TILE_W-STREL_SIZE+1)) &amp;&amp; (threadIdx.y &gt; (TILE_H-STREL_SIZE+1)))
             data[threadIdx.x + STREL_SIZE-1][threadIdx.y + STREL_SIZE-1] = buffer_in[(row + STREL_SIZE-1)*w + col + STREL_SIZE-1];
  
      //wait for all threads to finish read
       __syncthreads();
  
  
  
        unsigned char min_value = 255;
        for(x=0;x&lt;STREL_SIZE;x++){
            for(y=0;y&lt;STREL_SIZE;y++){
                min_value = min( (data[threadIdx.x+x][threadIdx.y+y]) , min_value);
                }
  
            }
        buffer_out[gLoc]= min_value;
        }
  
      }
  </code></pre>
  
  <p>my output now is:</p>
  
  <p><a href="https://i.stack.imgur.com/NUQhX.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/NUQhX.jpg" alt="new output"></a></p>
  
  <p><strong>UPDATED 2(version 2 -working-):</strong></p>
  
  <p>i have implemented another version of algorithm.To do that i follow that <a href="https://www.bsc.es/sites/default/files/public/computer_science/extreme_computing/bsc_intro_cuda_programming.pdf" rel="nofollow noreferrer">slide</a> that i found very useful and well explained,in particular the part in wich the author talk about convolution <strong>slide 27</strong>.</p>
  
  <p>i change the <strong>block</strong> and <strong>grid</strong> settings to :</p>
  
  <pre><code>dim3 block(20,20);
  dim3 grid(512/(block.x)+ block.x,512/(block.y)+block.y);
  </code></pre>
  
  <p>the kernel call instead ramain the same:</p>
  
  <pre><code>erode_multiple_img_SM&lt;&lt;&lt;grid,block&gt;&gt;&gt;(dimage_src,dimage_dst,512,512);
  </code></pre>
  
  <p>where the argument of the kernel are:</p>
  
  <ol>
  <li><strong>dimage_src:</strong> buffer of unsigned char with size <strong>height x width</strong> that contain input image.</li>
  <li><strong>dimage_dst:**buffer of unsigned char with size **height x width</strong> that contain output image, that my kernel produced.</li>
  <li><strong>512:</strong> the third argument is the width of the image.</li>
  <li><strong>512:</strong> the fourth argument is the height of the image.</li>
  </ol>
  
  <p><strong>remember</strong> my image sample are <strong>black and white</strong> but this version of erosion can work with <strong>grayscale</strong> too.</p>
  
  <p>here my <strong>working Kernel</strong>:</p>
  
  <pre><code>#define STREL_W 5
  #define STREL_H 5
  
  #define STREL_SIZE 5
  
  
  #define TILE_W 16
  #define TILE_H 16
  
  #define R (STREL_SIZE/2)
  
  
  #define BLOCK_W (TILE_W+(2*R))
  #define BLOCK_H (TILE_H+(2*R))
  
  __global__ void erode_multiple_img_working(unsigned char * buffer_in,
                              unsigned char * buffer_out,
                              int w,int h ){
  
  
      __shared__ unsigned char fast_acc_mat[BLOCK_w][BLOCK_H];
  
      int ty = threadIdx.y;
      int tx = threadIdx.x;
  
  
      int row_o = blockIdx.y * TILE_W + ty;
      int col_o = blockIdx.x * TILE_H + tx;
  
  
      int row_i = row_o - R;
      int col_i = col_o - R;
  
      //in of img size
      if((row_i &gt;= 0) &amp;&amp; (row_i &lt; h) &amp;&amp; (col_i &gt;= 0) &amp;&amp; (col_i &lt; w) ){
  
          fast_acc_mat[ty][tx] = buffer_in[ row_i * w + col_i];
  
      }
      else{
  
          fast_acc_mat[ty][tx] = 0;
  
      }
  
  
      __syncthreads();
  
  
  
  
  
      if( ty &lt; TILE_H &amp;&amp; tx &lt; TILE_W ){
  
          unsigned char min_val=255;
          for(int i = 0; i &lt; STREL_SIZE; i++) {
              for(int j = 0; j &lt; STREL_SIZE; j++) {
  
                  min_val = min( fast_acc_mat[i+ty][j+tx] , min_val );
  
              }
          }
          if(row_o &lt; h &amp;&amp; col_o &lt; w)
                  buffer_out[row_o * w + col_o] = min_val;
  
          }
  
       }
  </code></pre>
  
  <p>and this is my eroded image(<strong>output</strong>):</p>
  
  <p><a href="https://i.stack.imgur.com/qOeQT.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/qOeQT.jpg" alt="image eroded"></a></p>
  
  <p>I realized a scheme that show how the part of the algorithm  described by <strong>Eric</strong> load pixel of a <strong>TILE</strong> in shared memory :</p>
  
  <p><a href="https://i.stack.imgur.com/zaRbR.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/zaRbR.jpg" alt="enter image description here"></a></p>