paper_fast_median/code/kernMedianSeparable.cu~

   1 __global__ void kernel_medianV_sh( short *output, int i_dim, int j_dim, int r)
   2 {
   3
   4   int idc, val, min, max, inf, egal, sup, mxinf, minsup, estim ;
   5
   6   //coordinates in the block
   7   int ib = threadIdx.y ;
   8   int jb = threadIdx.x ;
   9   int idx_h = __mul24(ib+r,blockDim.x) + jb ;   // index pixel deans shmem (bloc+halo)
  10   int offset = __mul24(blockDim.x,r) ;
  11
  12   // coordonnees absolues du point
  13   int j = __mul24(blockIdx.x,blockDim.x) + jb ;
  14   int i = __mul24(blockIdx.y,blockDim.y) + ib ;
  15
  16   extern __shared__ int buff[] ;
  17   /***********************************************************************************
  18    *              CHARGEMENT DATA EN SHARED MEM
  19    ***********************************************************************************/
  20   buff[ idx_h ] = tex2D(tex_img_ins, j, i) ;
  21
  22   if (ib < r)
  23         {
  24           buff[ idx_h - offset ]    = tex2D(tex_img_ins, j, i-r) ;
  25         } else
  26         if (ib >= (blockDim.y-r))
  27           {
  28                 buff[ idx_h + offset  ]    = tex2D(tex_img_ins, j, i+r) ;
  29           }
  30
  31   __syncthreads() ;
  32   /**********************************************************************************************
  33    *               TRI VERTICAL par algo TORBEN MOGENSEN
  34    *          (a little bit slow but saves memory => faster !)
  35    **********************************************************************************************/
  36   min = max = buff[ ib*blockDim.x +jb] ;
  37
  38   for (idc= 0 ; idc< 2*r+1 ; idc++ )
  39         {
  40           val = buff[ __mul24(ib+idc, blockDim.x) +jb ] ;
  41           if ( val  < min ) min = val ;
  42           if ( val  > max ) max = val ;
  43         }
  44
  45   while (1)
  46         {
  47           estim = (min+max)/2 ;
  48           inf = sup = egal = 0  ;
  49           mxinf = min ;
  50           minsup= max ;
  51           for (idc =0; idc< 2*r+1 ; idc++)
  52                 {
  53                   val = buff[ __mul24(ib+idc, blockDim.x) +jb ] ;
  54                   if( val < estim )
  55                         {
  56                           inf++;
  57                           if( val > mxinf) mxinf = val ;
  58                         } else if (val > estim)
  59                         {
  60                           sup++;
  61                           if( val < minsup) minsup = val ;
  62                         } else egal++ ;
  63                 }
  64           if ( (inf <= (r+1))&&(sup <=(r+1)) ) break ;
  65           else if (inf>sup) max = mxinf ;
  66           else min = minsup ;
  67           }
  68
  69   if ( inf >= r+1 ) val = mxinf ;
  70   else if (inf+egal >= r+1) val = estim ;
  71   else val = minsup ;
  72
  73   output[ __mul24(j, i_dim)  +i  ] =  val ;
  74 }