Advertisement
Bisqwit

blur.hh

Sep 21st, 2017
262
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 2.96 KB | None | 0 0
  1. #include <cmath>
  2.  
  3. /* blur(): Really fast O(n) gaussian blur algorithm (gaussBlur_4)
  4.  * By Ivan Kuckir with ideas from Wojciech Jarosz
  5.  * Adapted from http://blog.ivank.net/fastest-gaussian-blur.html
  6.  *
  7.  * input:  The two-dimensional array of input signal. Must contain w*h elements.
  8.  * output: Where the two-dimensional array of blurred signal will be written
  9.  * temp:   Another array, for temporary use. Same size as input and output.
  10.  * w:      Width of array
  11.  * h:      Height of array.
  12.  * sigma:  Blurring kernel size. Must be smaller than w and h.
  13.  * n_boxes: Controls the blurring quality. 1 = box filter. 3 = pretty good filter.
  14.  *          Higher number = diminishingly better results, but linearly slower.
  15.  * elem_t: Type of elements. Should be integer type.
  16.  */
  17. template<unsigned n_boxes, typename elem_t>
  18. void blur(const elem_t* input, elem_t* output, elem_t* temp,
  19.           unsigned w,unsigned h,float sigma)
  20. {
  21.     auto wIdeal = std::sqrt((12*sigma*sigma/n_boxes)+1);  // Ideal averaging filter width
  22.     unsigned wl = wIdeal; if(wl%2==0) --wl;
  23.     unsigned wu = wl+2;
  24.     auto mIdeal = (12*sigma*sigma - n_boxes*wl*wl - 4*n_boxes*wl - 3*n_boxes)/(-4.*wl - 4);
  25.     unsigned m = std::round(mIdeal);
  26.     const elem_t* data = input;
  27.     for(unsigned n=0; n<n_boxes; ++n)
  28.     {
  29.         unsigned r = ((n<m ? wl : wu) - 1)/2; // IDK should this be float?
  30.         // boxBlur_4:
  31.         float iarr = 1.f / (r+r+1);
  32.         // boxBlurH_4 (blur horizontally for each row):
  33.         const elem_t* scl = data; elem_t* tcl = temp;
  34.         for(unsigned i=0; i<h; ++i)
  35.         {
  36.             auto ti = i*w, li = ti, ri = ti+r;
  37.             auto fv = scl[ti], lv = scl[ti+w-1]; int val = 0;
  38.             #pragma omp simd reduction(+:val)
  39.             for(unsigned j=0; j<r; j++) val += scl[ti+j];
  40.             val += (r+1)*fv;
  41.             for(unsigned j=0  ; j<=r ; j++) { val += scl[ri++] - fv       ;   tcl[ti++] = std::round(val*iarr); }
  42.             for(unsigned j=r+1; j<w-r; j++) { val += scl[ri++] - scl[li++];   tcl[ti++] = std::round(val*iarr); }
  43.             for(unsigned j=w-r; j<w  ; j++) { val += lv        - scl[li++];   tcl[ti++] = std::round(val*iarr); }
  44.         }
  45.         // boxBlurT_4 (blur vertically for each column)
  46.         scl = temp; tcl = output;
  47.         for(unsigned i=0; i<w; ++i)
  48.         {
  49.             auto ti = i, li = ti, ri = ti+r*w;
  50.             auto fv = scl[ti], lv = scl[ti+w*(h-1)]; int val = 0;
  51.             #pragma omp simd reduction(+:val)
  52.             for(unsigned j=0; j<r;  ++j) val += scl[ti + j*w];
  53.             val += (r+1)*fv;
  54.             for(unsigned j=0; j<=r; ++j)    { val += scl[ri] - fv     ;  tcl[ti] = std::round(val*iarr);  ri+=w; ti+=w; }
  55.             for(unsigned j=r+1; j<h-r; ++j) { val += scl[ri] - scl[li];  tcl[ti] = std::round(val*iarr);  li+=w; ri+=w; ti+=w; }
  56.             for(unsigned j=h-r; j<h; ++j)   { val += lv      - scl[li];  tcl[ti] = std::round(val*iarr);  li+=w; ti+=w; }
  57.         }
  58.         data = output;
  59.     }
  60. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement