include/cpuground.h

Go to the documentation of this file.

#ifndef CPUFILTER_H
#define CPUFILTER_H

//== INCLUDES =================================================================

#include <cmath>

#include <defs.h>
#include <extension.h>

//== NAMESPACES ===============================================================

namespace gpufilter {

//== IMPLEMENTATION ===========================================================

//-- Image --------------------------------------------------------------------

template< class T >
void extend_image( T *& ext_img,
                   int& ext_w,
                   int& ext_h,
                   const T *img,
                   const int& w,
                   const int& h,
                   const int& extb,
                   const initcond& ic ) {
    int bleft, btop, bright, bbottom;
    calc_borders( bleft, btop, bright, bbottom, w, h, extb );
    ext_w = w+bleft+bright;
    ext_h = h+btop+bbottom;
    ext_img = new float[ ext_w * ext_h ];
    for (int i = -btop; i < h+bbottom; ++i) {
        for (int j = -bleft; j < w+bright; ++j) {
            ext_img[(i+btop)*ext_w+(j+bleft)] = lookat(img, i, j, h, w, ic);
        }
    }
}

template< class T >
void extract_image( T *img,
                    const int& w,
                    const int& h,
                    T *& ext_img,
                    const int& ext_w,
                    const int& extb ) {
    int bleft, btop, bright, bbottom;
    calc_borders( bleft, btop, bright, bbottom, w, h, extb );
    for (int i = 0; i < h; ++i) {
        for (int j = 0; j < w; ++j) {
            img[i*w+j] = ext_img[(i+btop)*ext_w+(j+bleft)];
        }
    }
    delete [] ext_img;
}

//-- First-order Filter -------------------------------------------------------

template< class T >
void rcfr( T *inout,
           const int& w,
           const int& h,
           const T& b0,
           const T& a1,
           const bool& ff = false ) {
    for (int j = 0; j < w; j++) {
        int i;
        // p=0 is not initial condition based, it is due
        // to the filter order outside image + extension
        T p = (T)0;
        T c;
        for (i = 0; i < h; i++) {
            c = inout[i*w+j];
            p = c*b0 - p*a1;
            inout[i*w+j] = p;
        }
        if (ff) continue;
        p = (T)0;
        for (i--; i >= 0; i--) {
            c = inout[i*w+j];
            p = c*b0 - p*a1;
            inout[i*w+j] = p;
        }
    }
}

template< class T >
void rrfr( T *inout,
           const int& w,
           const int& h,
           const T& b0,
           const T& a1,
           const bool& ff = false ) {
    for (int i = 0; i < h; i++, inout += w) {
        int j;
        // p=0 is not initial condition based, it is due
        // to the filter order outside image + extension
        T p = (T)0;
        T c;
        for (j = 0; j < w; j++) {
            c = inout[j];
            p = c*b0 - p*a1;
            inout[j] = p;
        }
        if (ff) continue;
        p = (T)0;
        for (j--; j >= 0; j--) {
            c = inout[j];
            p = c*b0 - p*a1;
            inout[j] = p;
        }
    }
}
template< class T >
void r( T *inout,
        const int& w,
        const int& h,
        const T& b0,
        const T& a1,
        const bool& ff = false,
        const int& extb = 0,
        const initcond& ic = zero ) {
    if (extend(w, h, extb)) {
        int ext_w, ext_h;
        float *ext_inout;
        extend_image(ext_inout, ext_w, ext_h, inout, w, h, extb, ic);
        rcfr(ext_inout, ext_w, ext_h, b0, a1, ff);
        rrfr(ext_inout, ext_w, ext_h, b0, a1, ff);
        extract_image(inout, w, h, ext_inout, ext_w, extb);
    } else {
        rcfr(inout, w, h, b0, a1, ff);
        rrfr(inout, w, h, b0, a1, ff);
    }
}

//-- Second-order Filter ------------------------------------------------------

template< class T >
void rcfr( T *inout,
           const int& w,
           const int& h,
           const T& b0,
           const T& a1,
           const T& a2,
           const bool& ff = false ) {
    for (int j = 0; j < w; j++) {
        int i;
        // pp=p=0 is not initial condition based, it is due
        // to the filter order outside image + extension
        T pp = (T)0;
        T p = (T)0;
        T c;
        for (i = 0; i < h; i++) {
            c = inout[i*w+j];
            c = c*b0 - p*a1 - pp*a2;
            pp = p;
            p = c;
            inout[i*w+j] = p;
        }
        if (ff) continue;
        pp = p = (T)0;
        for (i--; i >= 0; i--) {
            c = inout[i*w+j];
            c = c*b0 - p*a1 - pp*a2;
            pp = p;
            p = c;
            inout[i*w+j] = p;
        }
    }
}

template< class T >
void rrfr( T *inout,
           const int& w,
           const int& h,
           const T& b0,
           const T& a1,
           const T& a2,
           const bool& ff = false ) {
    for (int i = 0; i < h; i++, inout += w) {
        int j;
        // pp=p=0 is not initial condition based, it is due
        // to the filter order outside image + extension
        T pp = (T)0;
        T p = (T)0;
        T c;
        for (j = 0; j < w; j++) {
            c = inout[j];
            c = c*b0 - p*a1 - pp*a2;
            pp = p;
            p = c;
            inout[j] = p;
        }
        if (ff) continue;
        pp = p = (T)0;
        for (j--; j >= 0; j--) {
            c = inout[j];
            c = c*b0 - p*a1 - pp*a2;
            pp = p;
            p = c;
            inout[j] = p;
        }
    }
}

template< class T >
void r( T *inout,
        const int& w,
        const int& h,
        const T& b0,
        const T& a1,
        const T& a2,
        const bool& ff = false,
        const int& extb = 0,
        const initcond& ic = zero ) {
    if (extend(w, h, extb)) {
        int ext_w, ext_h;
        float *ext_inout;
        extend_image(ext_inout, ext_w, ext_h, inout, w, h, extb, ic);
        rcfr(ext_inout, ext_w, ext_h, b0, a1, a2, ff);
        rrfr(ext_inout, ext_w, ext_h, b0, a1, a2, ff);
        extract_image(inout, w, h, ext_inout, ext_w, extb);
    } else {
        rcfr(inout, w, h, b0, a1, a2, ff);
        rrfr(inout, w, h, b0, a1, a2, ff);
    }
}

//-- SAT ----------------------------------------------------------------------

template< class T >
void sat_cpu( T *in,
              const int& w,
              const int& h ) {
    r(in, w, h, (T)1, (T)-1, true);
}
//-- Gaussian -----------------------------------------------------------------

template< class T >
void gaussian_cpu( T **in,
                   const int& w,
                   const int& h,
                   const int& depth,
                   const T& s,
                   const int& extb = 1,
                   const initcond& ic = clamp ) {
    T b10, a11, b20, a21, a22;
    weights1(s, b10, a11);
    weights2(s, b20, a21, a22);
    for (int c = 0; c < depth; c++) {
        r(in[c], w, h, b10, a11, false, extb, ic);
        r(in[c], w, h, b20, a21, a22, false, extb, ic);
    }
}
template< class T >
void gaussian_cpu( T *in,
                   const int& w,
                   const int& h,
                   const T& s,
                   const int& extb = 1,
                   const initcond& ic = clamp ) {
    T b10, a11, b20, a21, a22;
    weights1(s, b10, a11);
    weights2(s, b20, a21, a22);
    r(in, w, h, b10, a11, false, extb, ic);
    r(in, w, h, b20, a21, a22, false, extb, ic);
}
//-- BSpline ------------------------------------------------------------------

template< class T >
void bspline3i_cpu( T **in,
                    const int& w,
                    const int& h,
                    const int& depth,
                    const int& extb = 1,
                    const initcond& ic = mirror ) {
    const T alpha = (T)2 - sqrt((T)3);
    for (int c = 0; c < depth; c++) {
        r(in[c], w, h, (T)1+alpha, alpha, false, extb, ic);
    }
}

template< class T >
void bspline3i_cpu( T *in,
                    const int& w,
                    const int& h,
                    const int& extb = 1,
                    const initcond& ic = mirror ) {
    const T alpha = (T)2 - sqrt((T)3);
    r(in, w, h, (T)1+alpha, alpha, false, extb, ic);
}
//=============================================================================
} // namespace gpufilter
//=============================================================================
#endif // CPUFILTER_H
//=============================================================================