modules/ocl/src/opencl/filtering_filter2D.cl

   1 /*M///////////////////////////////////////////////////////////////////////////////////////
   2 //
   3 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
   4 //
   5 //  By downloading, copying, installing or using the software you agree to this license.
   6 //  If you do not agree to this license, do not download, install,
   7 //  copy or use the software.
   8 //
   9 //
  10 //                           License Agreement
  11 //                For Open Source Computer Vision Library
  12 //
  13 // Copyright (C) 2010-2013, Advanced Micro Devices, Inc., all rights reserved.
  14 // Third party copyrights are property of their respective owners.
  15 //
  16 // Redistribution and use in source and binary forms, with or without modification,
  17 // are permitted provided that the following conditions are met:
  18 //
  19 //   * Redistribution's of source code must retain the above copyright notice,
  20 //     this list of conditions and the following disclaimer.
  21 //
  22 //   * Redistribution's in binary form must reproduce the above copyright notice,
  23 //     this list of conditions and the following disclaimer in the documentation
  24 //     and/or other materials provided with the distribution.
  25 //
  26 //   * The name of the copyright holders may not be used to endorse or promote products
  27 //     derived from this software without specific prior written permission.
  28 //
  29 // This software is provided by the copyright holders and contributors as is and
  30 // any express or implied warranties, including, but not limited to, the implied
  31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
  32 // In no event shall the Intel Corporation or contributors be liable for any direct,
  33 // indirect, incidental, special, exemplary, or consequential damages
  34 // (including, but not limited to, procurement of substitute goods or services;
  35 // loss of use, data, or profits; or business interruption) however caused
  36 // and on any theory of liability, whether in contract, strict liability,
  37 // or tort (including negligence or otherwise) arising in any way out of
  38 // the use of this software, even if advised of the possibility of such damage.
  39 //
  40 //M*/
  41
  42 #ifdef BORDER_REPLICATE
  43 //BORDER_REPLICATE:     aaaaaa|abcdefgh|hhhhhhh
  44 #define ADDR_L(i, l_edge, r_edge)  ((i) <  (l_edge) ? (l_edge)   : (i))
  45 #define ADDR_R(i, r_edge, addr)    ((i) >= (r_edge) ? (r_edge)-1 : (addr))
  46 #define ADDR_H(i, t_edge, b_edge)  ((i) <  (t_edge) ? (t_edge)   :(i))
  47 #define ADDR_B(i, b_edge, addr)    ((i) >= (b_edge) ? (b_edge)-1 :(addr))
  48 #endif
  49
  50 #ifdef BORDER_REFLECT
  51 //BORDER_REFLECT:       fedcba|abcdefgh|hgfedcb
  52 #define ADDR_L(i, l_edge, r_edge)  ((i) <  (l_edge) ? -(i)-1               : (i))
  53 #define ADDR_R(i, r_edge, addr)    ((i) >= (r_edge) ? -(i)-1+((r_edge)<<1) : (addr))
  54 #define ADDR_H(i, t_edge, b_edge)  ((i) <  (t_edge) ? -(i)-1 : (i))
  55 #define ADDR_B(i, b_edge, addr)    ((i) >= (b_edge) ? -(i)-1+((b_edge)<<1) : (addr))
  56 #endif
  57
  58 #ifdef BORDER_REFLECT_101
  59 //BORDER_REFLECT_101:   gfedcb|abcdefgh|gfedcba
  60 #define ADDR_L(i, l_edge, r_edge)  ((i) <  (l_edge) ? -(i)                 : (i))
  61 #define ADDR_R(i, r_edge, addr)    ((i) >= (r_edge) ? -(i)-2+((r_edge)<<1) : (addr))
  62 #define ADDR_H(i, t_edge, b_edge)  ((i) <  (t_edge) ? -(i)                 : (i))
  63 #define ADDR_B(i, b_edge, addr)    ((i) >= (b_edge) ? -(i)-2+((b_edge)<<1) : (addr))
  64 #endif
  65
  66 //blur function does not support BORDER_WRAP
  67 #ifdef BORDER_WRAP
  68 //BORDER_WRAP:          cdefgh|abcdefgh|abcdefg
  69 #define ADDR_L(i, l_edge, r_edge)  ((i) <  (l_edge) ? (i)+(r_edge) : (i))
  70 #define ADDR_R(i, r_edge, addr)    ((i) >= (r_edge) ? (i)-(r_edge) : (addr))
  71 #define ADDR_H(i, t_edge, b_edge)  ((i) <  (t_edge) ? (i)+(b_edge) : (i))
  72 #define ADDR_B(i, b_edge, addr)    ((i) >= (b_edge) ? (i)-(b_edge) : (addr))
  73 #endif
  74
  75 #ifdef EXTRA_EXTRAPOLATION // border > src image size
  76 #ifdef BORDER_CONSTANT
  77 // None
  78 #elif defined BORDER_REPLICATE
  79 #define EXTRAPOLATE(x, y, minX, minY, maxX, maxY) \
  80     { \
  81         x = max(min(x, maxX - 1), minX); \
  82         y = max(min(y, maxY - 1), minY); \
  83     }
  84 #elif defined BORDER_WRAP
  85 #define EXTRAPOLATE(x, y, minX, minY, maxX, maxY) \
  86     { \
  87         if (x < minX) \
  88             x -= ((x - maxX + 1) / maxX) * maxX; \
  89         if (x >= maxX) \
  90             x %= maxX; \
  91         if (y < minY) \
  92             y -= ((y - maxY + 1) / maxY) * maxY; \
  93         if (y >= maxY) \
  94             y %= maxY; \
  95     }
  96 #elif defined(BORDER_REFLECT) || defined(BORDER_REFLECT_101)
  97 #define EXTRAPOLATE_(x, y, minX, minY, maxX, maxY, delta) \
  98     { \
  99         if (maxX - minX == 1) \
 100             x = minX; \
 101         else \
 102             do \
 103             { \
 104                 if (x < minX) \
 105                     x = -(x - minX) - 1 + delta; \
 106                 else \
 107                     x = maxX - 1 - (x - maxX) - delta; \
 108             } \
 109             while (x >= maxX || x < minX); \
 110         \
 111         if (maxY - minY == 1) \
 112             y = minY; \
 113         else \
 114             do \
 115             { \
 116                 if (y < minY) \
 117                     y = -(y - minY) - 1 + delta; \
 118                 else \
 119                     y = maxY - 1 - (y - maxY) - delta; \
 120             } \
 121             while (y >= maxY || y < minY); \
 122     }
 123 #ifdef BORDER_REFLECT
 124 #define EXTRAPOLATE(x, y, minX, minY, maxX, maxY) EXTRAPOLATE_(x, y, minX, minY, maxX, maxY, 0)
 125 #elif defined(BORDER_REFLECT_101)
 126 #define EXTRAPOLATE(x, y, minX, minY, maxX, maxY) EXTRAPOLATE_(x, y, minX, minY, maxX, maxY, 1)
 127 #endif
 128 #else
 129 #error No extrapolation method
 130 #endif
 131 #else
 132 #define EXTRAPOLATE(x, y, minX, minY, maxX, maxY) \
 133     { \
 134         int _row = y - minY, _col = x - minX; \
 135         _row = ADDR_H(_row, 0, maxY - minY); \
 136         _row = ADDR_B(_row, maxY - minY, _row); \
 137         y = _row + minY; \
 138         \
 139         _col = ADDR_L(_col, 0, maxX - minX); \
 140         _col = ADDR_R(_col, maxX - minX, _col); \
 141         x = _col + minX; \
 142     }
 143 #endif
 144
 145 #if USE_DOUBLE
 146 #ifdef cl_amd_fp64
 147 #pragma OPENCL EXTENSION cl_amd_fp64:enable
 148 #elif defined (cl_khr_fp64)
 149 #pragma OPENCL EXTENSION cl_khr_fp64:enable
 150 #endif
 151 #define FPTYPE double
 152 #define CONVERT_TO_FPTYPE CAT(convert_double, VEC_SIZE)
 153 #else
 154 #define FPTYPE float
 155 #define CONVERT_TO_FPTYPE CAT(convert_float, VEC_SIZE)
 156 #endif
 157
 158 #if DATA_DEPTH == 0
 159 #define BASE_TYPE uchar
 160 #elif DATA_DEPTH == 1
 161 #define BASE_TYPE char
 162 #elif DATA_DEPTH == 2
 163 #define BASE_TYPE ushort
 164 #elif DATA_DEPTH == 3
 165 #define BASE_TYPE short
 166 #elif DATA_DEPTH == 4
 167 #define BASE_TYPE int
 168 #elif DATA_DEPTH == 5
 169 #define BASE_TYPE float
 170 #elif DATA_DEPTH == 6
 171 #define BASE_TYPE double
 172 #else
 173 #error data_depth
 174 #endif
 175
 176 #define __CAT(x, y) x##y
 177 #define CAT(x, y) __CAT(x, y)
 178
 179 #define uchar1 uchar
 180 #define char1 char
 181 #define ushort1 ushort
 182 #define short1 short
 183 #define int1 int
 184 #define float1 float
 185 #define double1 double
 186
 187 #define convert_uchar1_sat_rte convert_uchar_sat_rte
 188 #define convert_char1_sat_rte convert_char_sat_rte
 189 #define convert_ushort1_sat_rte convert_ushort_sat_rte
 190 #define convert_short1_sat_rte convert_short_sat_rte
 191 #define convert_int1_sat_rte convert_int_sat_rte
 192 #define convert_float1
 193 #define convert_double1
 194
 195 #if DATA_DEPTH == 5 || DATA_DEPTH == 6
 196 #define CONVERT_TO_TYPE CAT(CAT(convert_, BASE_TYPE), VEC_SIZE)
 197 #else
 198 #define CONVERT_TO_TYPE CAT(CAT(CAT(convert_, BASE_TYPE), VEC_SIZE), _sat_rte)
 199 #endif
 200
 201 #define VEC_SIZE DATA_CHAN
 202
 203 #define VEC_TYPE CAT(BASE_TYPE, VEC_SIZE)
 204 #define TYPE VEC_TYPE
 205
 206 #define SCALAR_TYPE CAT(FPTYPE, VEC_SIZE)
 207
 208 #define INTERMEDIATE_TYPE CAT(FPTYPE, VEC_SIZE)
 209
 210 struct RectCoords
 211 {
 212     int x1, y1, x2, y2;
 213 };
 214
 215 //#define DEBUG
 216 #ifdef DEBUG
 217 #define DEBUG_ONLY(x) x
 218 #define ASSERT(condition) do { if (!(condition)) { printf("BUG in boxFilter kernel (global=%d,%d): " #condition "\n", get_global_id(0), get_global_id(1)); } } while (0)
 219 #else
 220 #define DEBUG_ONLY(x) (void)0
 221 #define ASSERT(condition) (void)0
 222 #endif
 223
 224
 225 inline INTERMEDIATE_TYPE readSrcPixel(int2 pos, __global TYPE *src, const unsigned int srcStepBytes, const struct RectCoords srcCoords
 226 #ifdef BORDER_CONSTANT
 227                , SCALAR_TYPE borderValue
 228 #endif
 229     )
 230 {
 231 #ifdef BORDER_ISOLATED
 232     if(pos.x >= srcCoords.x1 && pos.y >= srcCoords.y1 && pos.x < srcCoords.x2 && pos.y < srcCoords.y2)
 233 #else
 234     if(pos.x >= 0 && pos.y >= 0 && pos.x < srcCoords.x2 && pos.y < srcCoords.y2)
 235 #endif
 236     {
 237         __global TYPE* ptr = (__global TYPE*)((__global char*)src + pos.x * sizeof(TYPE) + pos.y * srcStepBytes);
 238         return CONVERT_TO_FPTYPE(*ptr);
 239     }
 240     else
 241     {
 242 #ifdef BORDER_CONSTANT
 243         return borderValue;
 244 #else
 245         int selected_col = pos.x;
 246         int selected_row = pos.y;
 247
 248         EXTRAPOLATE(selected_col, selected_row,
 249 #ifdef BORDER_ISOLATED
 250                 srcCoords.x1, srcCoords.y1,
 251 #else
 252                 0, 0,
 253 #endif
 254                 srcCoords.x2, srcCoords.y2
 255          );
 256
 257         // debug border mapping
 258         //printf("pos=%d,%d --> %d, %d\n", pos.x, pos.y, selected_col, selected_row);
 259
 260         pos = (int2)(selected_col, selected_row);
 261         if(pos.x >= 0 && pos.y >= 0 && pos.x < srcCoords.x2 && pos.y < srcCoords.y2)
 262         {
 263             __global TYPE* ptr = (__global TYPE*)((__global char*)src + pos.x * sizeof(TYPE) + pos.y * srcStepBytes);
 264             return CONVERT_TO_FPTYPE(*ptr);
 265         }
 266         else
 267         {
 268             // for debug only
 269             DEBUG_ONLY(printf("BUG in boxFilter kernel\n"));
 270             return (FPTYPE)(0.0f);
 271         }
 272 #endif
 273     }
 274 }
 275
 276 // INPUT PARAMETER: BLOCK_SIZE_Y (via defines)
 277
 278 __kernel
 279 __attribute__((reqd_work_group_size(LOCAL_SIZE, 1, 1)))
 280 void filter2D(__global TYPE *src, const unsigned int srcStepBytes, const int4 srcRC,
 281               __global TYPE *dst, const unsigned int dstStepBytes, const int4 dstRC,
 282 #ifdef BORDER_CONSTANT
 283               SCALAR_TYPE borderValue,
 284 #endif
 285               __constant FPTYPE* kernelData // transposed: [KERNEL_SIZE_X][KERNEL_SIZE_Y2_ALIGNED]
 286               )
 287 {
 288     const struct RectCoords srcCoords = {srcRC.s0, srcRC.s1, srcRC.s2, srcRC.s3}; // for non-isolated border: offsetX, offsetY, wholeX, wholeY
 289     struct RectCoords dstCoords = {dstRC.s0, dstRC.s1, dstRC.s2, dstRC.s3};
 290
 291     const int local_id = get_local_id(0);
 292     const int x = local_id + (LOCAL_SIZE - (KERNEL_SIZE_X - 1)) * get_group_id(0) - ANCHOR_X;
 293     const int y = get_global_id(1) * BLOCK_SIZE_Y;
 294
 295     INTERMEDIATE_TYPE data[KERNEL_SIZE_Y];
 296     __local INTERMEDIATE_TYPE sumOfCols[LOCAL_SIZE];
 297
 298     int2 srcPos = (int2)(srcCoords.x1 + x, srcCoords.y1 + y - ANCHOR_Y);
 299
 300     int2 pos = (int2)(dstCoords.x1 + x, dstCoords.y1 + y);
 301     __global TYPE* dstPtr = (__global TYPE*)((__global char*)dst + pos.x * sizeof(TYPE) + pos.y * dstStepBytes); // Pointer can be out of bounds!
 302     bool writeResult = (local_id >= ANCHOR_X && local_id < LOCAL_SIZE - (KERNEL_SIZE_X - 1 - ANCHOR_X) &&
 303                         pos.x >= dstCoords.x1 && pos.x < dstCoords.x2);
 304
 305 #if BLOCK_SIZE_Y > 1
 306     bool readAllpixels = true;
 307     int sy_index = 0; // current index in data[] array
 308
 309     dstCoords.y2 = min(dstCoords.y2, pos.y + BLOCK_SIZE_Y);
 310     for (;
 311          pos.y < dstCoords.y2;
 312          pos.y++,
 313          dstPtr = (__global TYPE*)((__global char*)dstPtr + dstStepBytes))
 314 #endif
 315     {
 316         ASSERT(pos.y < dstCoords.y2);
 317
 318         for (
 319 #if BLOCK_SIZE_Y > 1
 320             int sy = readAllpixels ? 0 : -1; sy < (readAllpixels ? KERNEL_SIZE_Y : 0);
 321 #else
 322             int sy = 0, sy_index = 0; sy < KERNEL_SIZE_Y;
 323 #endif
 324             sy++, srcPos.y++)
 325         {
 326             data[sy + sy_index] = readSrcPixel(srcPos, src, srcStepBytes, srcCoords
 327 #ifdef BORDER_CONSTANT
 328                     , borderValue
 329 #endif
 330                     );
 331         }
 332
 333         INTERMEDIATE_TYPE total_sum = 0;
 334         for (int sx = 0; sx < KERNEL_SIZE_X; sx++)
 335         {
 336             {
 337                 __constant FPTYPE* k = &kernelData[KERNEL_SIZE_Y2_ALIGNED * sx
 338 #if BLOCK_SIZE_Y > 1
 339                                                    + KERNEL_SIZE_Y - sy_index
 340 #endif
 341                                                    ];
 342                 INTERMEDIATE_TYPE tmp_sum = 0;
 343                 for (int sy = 0; sy < KERNEL_SIZE_Y; sy++)
 344                 {
 345                     tmp_sum += data[sy] * k[sy];
 346                 }
 347
 348                 sumOfCols[local_id] = tmp_sum;
 349                 barrier(CLK_LOCAL_MEM_FENCE);
 350             }
 351
 352             int id = local_id + sx - ANCHOR_X;
 353             if (id >= 0 && id < LOCAL_SIZE)
 354                total_sum += sumOfCols[id];
 355
 356             barrier(CLK_LOCAL_MEM_FENCE);
 357         }
 358
 359         if (writeResult)
 360         {
 361             ASSERT(pos.y >= dstCoords.y1 && pos.y < dstCoords.y2);
 362             *dstPtr = CONVERT_TO_TYPE(total_sum);
 363         }
 364
 365 #if BLOCK_SIZE_Y > 1
 366         readAllpixels = false;
 367 #if BLOCK_SIZE_Y > KERNEL_SIZE_Y
 368         sy_index = (sy_index + 1 <= KERNEL_SIZE_Y) ? sy_index + 1 : 1;
 369 #else
 370         sy_index++;
 371 #endif
 372 #endif // BLOCK_SIZE_Y == 1
 373     }
 374 }