快速双边滤波 附完整C代码

2017-10-05 21:25:23来源:cnblogs.com作者:落羽の殇人点击

分享

很早之前写过《双边滤波算法的简易实现bilateralFilter》。

当时学习参考的代码来自cuda的样例。

相关代码可以参阅:

https://github.com/johng12/cudaSamples/tree/master/cudaSamples/3_Imaging/bilateralFilter

由于算法逻辑非常清晰,就不多解释了。

需要补课的,请移步《o(1)复杂度之双边滤波算法的原理、流程、实现及效果。》

代码见:bilateralFilter_cpu.cpp 文件。

#include <math.h>#include <string.h>////////////////////////////////////////////////////////////////////////////////// export C interface#define EPSILON 1e-3extern "C" void updateGaussianGold(float delta, int radius);extern "C" void bilateralFilterGold(unsigned int *pSrc,                                    unsigned int *pDest,                                    float e_d,                                    int w, int h, int r);//variablesfloat gaussian[50];struct float4{    float x;    float y;    float z;    float w;    float4() {};    float4(float value)    {        x = y = z = w = value;    }};void updateGaussianGold(float delta, int radius){    for (int i = 0; i < 2 * radius + 1; i++)    {        int x = i - radius;        gaussian[i] = exp(-(x * x) /                          (2 * delta * delta));    }}float heuclideanLen(float4 a, float4 b, float d){    float mod = (b.x - a.x) * (b.x - a.x) +                (b.y - a.y) * (b.y - a.y) +                (b.z - a.z) * (b.z - a.z) +                (b.w - a.w) * (b.w - a.w);    return exp(-mod / (2 * d * d));}unsigned int hrgbaFloatToInt(float4 rgba){    unsigned int w = (((unsigned int)(fabs(rgba.w) * 255.0f)) & 0xff) << 24;    unsigned int z = (((unsigned int)(fabs(rgba.z) * 255.0f)) & 0xff) << 16;    unsigned int y = (((unsigned int)(fabs(rgba.y) * 255.0f)) & 0xff) << 8;    unsigned int x = ((unsigned int)(fabs(rgba.x) * 255.0f)) & 0xff;    return (w | z | y | x);}float4 hrgbaIntToFloat(unsigned int c){    float4 rgba;    rgba.x = (c & 0xff) * 0.003921568627f;       //  /255.0f;    rgba.y = ((c>>8) & 0xff) * 0.003921568627f;  //  /255.0f;    rgba.z = ((c>>16) & 0xff) * 0.003921568627f; //  /255.0f;    rgba.w = ((c>>24) & 0xff) * 0.003921568627f; //  /255.0f;    return rgba;}float4 mul(float a, float4 b){    float4 ans;    ans.x = a * b.x;    ans.y = a * b.y;    ans.z = a * b.z;    ans.w = a * b.w;    return ans;}float4 add4(float4 a, float4 b){    float4 ans;    ans.x = a.x + b.x;    ans.y = a.y + b.y;    ans.z = a.z + b.z;    ans.w = a.w + b.w;    return ans;}void bilateralFilterGold(unsigned int *pSrc,                         unsigned int *pDest,                         float e_d,                         int w, int h, int r){    float4 *hImage = new float4[w * h];    float domainDist, colorDist, factor;    for (int y = 0; y < h; y++)    {        for (int x = 0; x < w; x++)        {            hImage[y * w + x] = hrgbaIntToFloat(pSrc[y * w + x]);        }    }    for (int y = 0; y < h; y++)    {        for (int x = 0; x < w; x++)        {            float4 t(0.0f);            float sum = 0.0f;            for (int i = -r; i <= r; i++)            {                int neighborY = y + i;                //clamp the neighbor pixel, prevent overflow                if (neighborY < 0)                {                    neighborY = 0;                }                else if (neighborY >= h)                {                    neighborY = h - 1;                }                for (int j = -r; j <= r; j++)                {                    domainDist = gaussian[r + i] * gaussian[r + j];                    //clamp the neighbor pixel, prevent overflow                    int neighborX = x + j;                    if (neighborX < 0)                    {                        neighborX = 0;                    }                    else if (neighborX >= w)                    {                        neighborX = w - 1;                    }                    colorDist = heuclideanLen(hImage[neighborY * w + neighborX], hImage[y * w + x], e_d);                    factor = domainDist * colorDist;                    sum += factor;                    t = add4(t, mul(factor, hImage[neighborY * w + neighborX]));                }            }            pDest[y * w + x] = hrgbaFloatToInt(mul(1 / sum, t));        }    }    delete[] hImage;}

在网上很多人都采用双边滤波作为磨皮算法,包括比较知名的gpuimage库。

用的也是小半径的双边,而双边滤波有很多快速算法的变种。

详情可参阅:http://people.csail.mit.edu/sparis/bf/

而就在几个月前,github上有个朋友放出来了一份快速双边滤波的tiny版本代码。

https://github.com/ufoym/RecursiveBF/

代码见:rbf.hpp

#ifndef INCLUDE_RBF#define INCLUDE_RBF#include <math.h>#include <string.h>#define QX_DEF_CHAR_MAX 255/* ======================================================================RecursiveBF: A lightweight library for recursive bilateral filtering.-------------------------------------------------------------------------Intro:      Recursive bilateral filtering (developed by Qingxiong Yang)             is pretty fast compared with most edge-preserving filtering             methods.            -   computational complexity is linear in both input size and                 dimensionality            -   takes about 43 ms to process a one mega-pixel color image                (i7 1.8GHz & 4GB memory)            -   about 18x faster than Fast high-dimensional filtering                 using the permutohedral lattice            -   about 86x faster than Gaussian kd-trees for fast high-                dimensional filteringUsage:      // ----------------------------------------------------------            // Basic Usage            // ----------------------------------------------------------            unsigned char * img = ...;                    // input image            unsigned char * img_out = 0;            // output image            int width = ..., height = ..., channel = ...; // image size            recursive_bf(img, img_out,                          sigma_spatial, sigma_range,                          width, height, channel);            // ----------------------------------------------------------            // Advanced: using external buffer for better performance            // ----------------------------------------------------------            unsigned char * img = ...;                    // input image            unsigned char * img_out = 0;            // output image            int width = ..., height = ..., channel = ...; // image size            float * buffer = new float[                   // external buf                                 ( width * height* channel                                  + width * height                                 + width * channel                                  + width) * 2];            recursive_bf(img, img_out,                          sigma_spatial, sigma_range,                          width, height, channel,                          buffer);            delete[] buffer;Notice:     Large sigma_spatial/sigma_range parameter may results in             visible artifact which can be removed by an additional             filter with small sigma_spatial/sigma_range parameter.-------------------------------------------------------------------------Reference:  Qingxiong Yang, Recursive Bilateral Filtering,            European Conference on Computer Vision (ECCV) 2012, 399-413.====================================================================== */inline void recursive_bf(    unsigned char * img_in,     unsigned char *& img_out,     float sigma_spatial, float sigma_range,     int width, int height, int channel,     float * buffer /*= 0*/);// ----------------------------------------------------------------------inline void _recursive_bf(    unsigned char * img,    float sigma_spatial, float sigma_range,     int width, int height, int channel,    float * buffer = 0){    const int width_height = width * height;    const int width_channel = width * channel;    const int width_height_channel = width * height * channel;    bool is_buffer_internal = (buffer == 0);    if (is_buffer_internal)        buffer = new float[(width_height_channel + width_height                             + width_channel + width) * 2];    float * img_out_f = buffer;    float * img_temp = &img_out_f[width_height_channel];    float * map_factor_a = &img_temp[width_height_channel];    float * map_factor_b = &map_factor_a[width_height];     float * slice_factor_a = &map_factor_b[width_height];    float * slice_factor_b = &slice_factor_a[width_channel];    float * line_factor_a = &slice_factor_b[width_channel];    float * line_factor_b = &line_factor_a[width];        //compute a lookup table    float range_table[QX_DEF_CHAR_MAX + 1];    float inv_sigma_range = 1.0f / (sigma_range * QX_DEF_CHAR_MAX);    for (int i = 0; i <= QX_DEF_CHAR_MAX; i++)         range_table[i] = static_cast<float>(exp(-i * inv_sigma_range));    float alpha = static_cast<float>(exp(-sqrt(2.0) / (sigma_spatial * width)));    float ypr, ypg, ypb, ycr, ycg, ycb;    float fp, fc;    float inv_alpha_ = 1 - alpha;    for (int y = 0; y < height; y++)    {        float * temp_x = &img_temp[y * width_channel];        unsigned char * in_x = &img[y * width_channel];        unsigned char * texture_x = &img[y * width_channel];        *temp_x++ = ypr = *in_x++;         *temp_x++ = ypg = *in_x++;         *temp_x++ = ypb = *in_x++;        unsigned char tpr = *texture_x++;         unsigned char tpg = *texture_x++;        unsigned char tpb = *texture_x++;        float * temp_factor_x = &map_factor_a[y * width];        *temp_factor_x++ = fp = 1;        // from left to right        for (int x = 1; x < width; x++)         {            unsigned char tcr = *texture_x++;             unsigned char tcg = *texture_x++;             unsigned char tcb = *texture_x++;            unsigned char dr = abs(tcr - tpr);            unsigned char dg = abs(tcg - tpg);            unsigned char db = abs(tcb - tpb);            int range_dist = (((dr << 1) + dg + db) >> 2);            float weight = range_table[range_dist];            float alpha_ = weight*alpha;            *temp_x++ = ycr = inv_alpha_*(*in_x++) + alpha_*ypr;             *temp_x++ = ycg = inv_alpha_*(*in_x++) + alpha_*ypg;             *temp_x++ = ycb = inv_alpha_*(*in_x++) + alpha_*ypb;            tpr = tcr; tpg = tcg; tpb = tcb;            ypr = ycr; ypg = ycg; ypb = ycb;            *temp_factor_x++ = fc = inv_alpha_ + alpha_*fp;            fp = fc;        }        *--temp_x; *temp_x = 0.5f*((*temp_x) + (*--in_x));        *--temp_x; *temp_x = 0.5f*((*temp_x) + (*--in_x));        *--temp_x; *temp_x = 0.5f*((*temp_x) + (*--in_x));        tpr = *--texture_x;         tpg = *--texture_x;         tpb = *--texture_x;        ypr = *in_x; ypg = *in_x; ypb = *in_x;        *--temp_factor_x; *temp_factor_x = 0.5f*((*temp_factor_x) + 1);        fp = 1;        // from right to left        for (int x = width - 2; x >= 0; x--)         {            unsigned char tcr = *--texture_x;             unsigned char tcg = *--texture_x;             unsigned char tcb = *--texture_x;            unsigned char dr = abs(tcr - tpr);            unsigned char dg = abs(tcg - tpg);            unsigned char db = abs(tcb - tpb);            int range_dist = (((dr << 1) + dg + db) >> 2);            float weight = range_table[range_dist];            float alpha_ = weight * alpha;            ycr = inv_alpha_ * (*--in_x) + alpha_ * ypr;             ycg = inv_alpha_ * (*--in_x) + alpha_ * ypg;             ycb = inv_alpha_ * (*--in_x) + alpha_ * ypb;            *--temp_x; *temp_x = 0.5f*((*temp_x) + ycr);            *--temp_x; *temp_x = 0.5f*((*temp_x) + ycg);            *--temp_x; *temp_x = 0.5f*((*temp_x) + ycb);            tpr = tcr; tpg = tcg; tpb = tcb;            ypr = ycr; ypg = ycg; ypb = ycb;            fc = inv_alpha_ + alpha_*fp;            *--temp_factor_x;             *temp_factor_x = 0.5f*((*temp_factor_x) + fc);            fp = fc;        }    }    alpha = static_cast<float>(exp(-sqrt(2.0) / (sigma_spatial * height)));    inv_alpha_ = 1 - alpha;    float * ycy, * ypy, * xcy;    unsigned char * tcy, * tpy;    memcpy(img_out_f, img_temp, sizeof(float)* width_channel);    float * in_factor = map_factor_a;    float*ycf, *ypf, *xcf;    memcpy(map_factor_b, in_factor, sizeof(float) * width);    for (int y = 1; y < height; y++)    {        tpy = &img[(y - 1) * width_channel];        tcy = &img[y * width_channel];        xcy = &img_temp[y * width_channel];        ypy = &img_out_f[(y - 1) * width_channel];        ycy = &img_out_f[y * width_channel];        xcf = &in_factor[y * width];        ypf = &map_factor_b[(y - 1) * width];        ycf = &map_factor_b[y * width];        for (int x = 0; x < width; x++)        {            unsigned char dr = abs((*tcy++) - (*tpy++));            unsigned char dg = abs((*tcy++) - (*tpy++));            unsigned char db = abs((*tcy++) - (*tpy++));            int range_dist = (((dr << 1) + dg + db) >> 2);            float weight = range_table[range_dist];            float alpha_ = weight*alpha;            for (int c = 0; c < channel; c++)                 *ycy++ = inv_alpha_*(*xcy++) + alpha_*(*ypy++);            *ycf++ = inv_alpha_*(*xcf++) + alpha_*(*ypf++);        }    }    int h1 = height - 1;    ycf = line_factor_a;    ypf = line_factor_b;    memcpy(ypf, &in_factor[h1 * width], sizeof(float) * width);    for (int x = 0; x < width; x++)         map_factor_b[h1 * width + x] = 0.5f*(map_factor_b[h1 * width + x] + ypf[x]);    ycy = slice_factor_a;    ypy = slice_factor_b;    memcpy(ypy, &img_temp[h1 * width_channel], sizeof(float)* width_channel);    int k = 0;     for (int x = 0; x < width; x++) {        for (int c = 0; c < channel; c++) {            int idx = (h1 * width + x) * channel + c;            img_out_f[idx] = 0.5f*(img_out_f[idx] + ypy[k++]) / map_factor_b[h1 * width + x];        }    }    for (int y = h1 - 1; y >= 0; y--)    {        tpy = &img[(y + 1) * width_channel];        tcy = &img[y * width_channel];        xcy = &img_temp[y * width_channel];        float*ycy_ = ycy;        float*ypy_ = ypy;        float*out_ = &img_out_f[y * width_channel];        xcf = &in_factor[y * width];        float*ycf_ = ycf;        float*ypf_ = ypf;        float*factor_ = &map_factor_b[y * width];        for (int x = 0; x < width; x++)        {            unsigned char dr = abs((*tcy++) - (*tpy++));            unsigned char dg = abs((*tcy++) - (*tpy++));            unsigned char db = abs((*tcy++) - (*tpy++));            int range_dist = (((dr << 1) + dg + db) >> 2);            float weight = range_table[range_dist];            float alpha_ = weight*alpha;            float fcc = inv_alpha_*(*xcf++) + alpha_*(*ypf_++);            *ycf_++ = fcc;            *factor_ = 0.5f * (*factor_ + fcc);            for (int c = 0; c < channel; c++)            {                float ycc = inv_alpha_*(*xcy++) + alpha_*(*ypy_++);                *ycy_++ = ycc;                *out_ = 0.5f * (*out_ + ycc) / (*factor_);                *out_++;            }            *factor_++;        }        memcpy(ypy, ycy, sizeof(float) * width_channel);        memcpy(ypf, ycf, sizeof(float) * width);    }    for (int i = 0; i < width_height_channel; ++i)        img[i] = static_cast<unsigned char>(img_out_f[i]);    if (is_buffer_internal)        delete[] buffer;}inline void recursive_bf(    unsigned char * img_in,    unsigned char *& img_out,    float sigma_spatial, float sigma_range,    int width, int height, int channel,    float * buffer = 0){    if (img_out == 0)        img_out = new unsigned char[width * height * channel];    for (int i = 0; i < width * height * channel; ++i)        img_out[i] = img_in[i];    _recursive_bf(img_out, sigma_spatial, sigma_range, width, height, channel, buffer);}#endif // INCLUDE_RBF

而这份代码的作者,真的很粗心,代码中有多处逻辑写错,代码风格也有待加强。

当然也有,有心者,在此基础上做了修正。

https://github.com/Fig1024/OP_RBF

这个作者还在这个基础上做了simd指令集优化。

优化代码暂时不看。

关注下基础实现的代码。即RBFilterPlain.cpp:

#include "stdafx.h"#include "RBFilterPlain.h"#include "stdafx.h"#include <algorithm>using namespace std;#define QX_DEF_CHAR_MAX 255CRBFilterPlain::CRBFilterPlain(){}CRBFilterPlain::~CRBFilterPlain(){    releaseMemory();}// assumes 3/4 channel images, 1 byte per channelvoid CRBFilterPlain::reserveMemory(int max_width, int max_height, int channels){    // basic sanity check    _ASSERT(max_width >= 10 && max_width < 10000);    _ASSERT(max_height >= 10 && max_height < 10000);    _ASSERT(channels >= 1 && channels <= 4);    releaseMemory();    m_reserve_width = max_width;    m_reserve_height = max_height;    m_reserve_channels = channels;    int width_height = m_reserve_width * m_reserve_height;    int width_height_channel = width_height * m_reserve_channels;    m_left_pass_color = new float[width_height_channel];    m_left_pass_factor = new float[width_height];    m_right_pass_color = new float[width_height_channel];    m_right_pass_factor = new float[width_height];    m_down_pass_color = new float[width_height_channel];    m_down_pass_factor = new float[width_height];    m_up_pass_color = new float[width_height_channel];    m_up_pass_factor = new float[width_height];}void CRBFilterPlain::releaseMemory(){    m_reserve_width = 0;    m_reserve_height = 0;    m_reserve_channels = 0;    if (m_left_pass_color)    {        delete[] m_left_pass_color;        m_left_pass_color = nullptr;    }    if (m_left_pass_factor)    {        delete[] m_left_pass_factor;        m_left_pass_factor = nullptr;    }    if (m_right_pass_color)    {        delete[] m_right_pass_color;        m_right_pass_color = nullptr;    }    if (m_right_pass_factor)    {        delete[] m_right_pass_factor;        m_right_pass_factor = nullptr;    }    if (m_down_pass_color)    {        delete[] m_down_pass_color;        m_down_pass_color = nullptr;    }    if (m_down_pass_factor)    {        delete[] m_down_pass_factor;        m_down_pass_factor = nullptr;    }    if (m_up_pass_color)    {        delete[] m_up_pass_color;        m_up_pass_color = nullptr;    }    if (m_up_pass_factor)    {        delete[] m_up_pass_factor;        m_up_pass_factor = nullptr;    }}int CRBFilterPlain::getDiffFactor(const unsigned char* color1, const unsigned char* color2) const{    int final_diff;    int component_diff[4];    // find absolute difference between each component    for (int i = 0; i < m_reserve_channels; i++)    {        component_diff[i] = abs(color1[i] - color2[i]);    }    // based on number of components, produce a single difference value in the 0-255 range    switch (m_reserve_channels)    {    case 1:        final_diff = component_diff[0];        break;    case 2:        final_diff = ((component_diff[0] + component_diff[1]) >> 1);        break;    case 3:        final_diff = ((component_diff[0] + component_diff[2]) >> 2) + (component_diff[1] >> 1);        break;    case 4:        final_diff = ((component_diff[0] + component_diff[1] + component_diff[2] + component_diff[3]) >> 2);        break;    default:        final_diff = 0;    }    _ASSERT(final_diff >= 0 && final_diff <= 255);    return final_diff;}// memory must be reserved before calling image filter// this implementation of filter uses plain C++, single threaded// channel count must be 3 or 4 (alpha not used)void CRBFilterPlain::filter(unsigned char* img_src, unsigned char* img_dst,    float sigma_spatial, float sigma_range,    int width, int height, int channel){    _ASSERT(img_src);    _ASSERT(img_dst);    _ASSERT(m_reserve_channels == channel);    _ASSERT(m_reserve_width >= width);    _ASSERT(m_reserve_height >= height);    // compute a lookup table    float alpha_f = static_cast<float>(exp(-sqrt(2.0) / (sigma_spatial * 255)));    float inv_alpha_f = 1.f - alpha_f;    float range_table_f[QX_DEF_CHAR_MAX + 1];    float inv_sigma_range = 1.0f / (sigma_range * QX_DEF_CHAR_MAX);    {        float ii = 0.f;        for (int i = 0; i <= QX_DEF_CHAR_MAX; i++, ii -= 1.f)        {            range_table_f[i] = alpha_f * exp(ii * inv_sigma_range);        }    }    ///////////////    // Left pass    {        const unsigned char* src_color = img_src;        float* left_pass_color = m_left_pass_color;        float* left_pass_factor = m_left_pass_factor;        for (int y = 0; y < height; y++)        {            const unsigned char* src_prev = src_color;            const float* prev_factor = left_pass_factor;            const float* prev_color = left_pass_color;            // process 1st pixel separately since it has no previous            *left_pass_factor++ = 1.f;            for (int c = 0; c < channel; c++)            {                *left_pass_color++ = *src_color++;            }            // handle other pixels            for (int x = 1; x < width; x++)            {                // determine difference in pixel color between current and previous                // calculation is different depending on number of channels                int diff = getDiffFactor(src_color, src_prev);                src_prev = src_color;                float alpha_f = range_table_f[diff];                *left_pass_factor++ = inv_alpha_f + alpha_f * (*prev_factor++);                for (int c = 0; c < channel; c++)                {                    *left_pass_color++ = inv_alpha_f * (*src_color++) + alpha_f * (*prev_color++);                }            }        }    }    ///////////////    // Right pass    {        // start from end and then go up to begining         int last_index = width * height * channel - 1;        const unsigned char* src_color = img_src + last_index;        float* right_pass_color = m_right_pass_color + last_index;        float* right_pass_factor = m_right_pass_factor + width * height - 1;        for (int y = 0; y < height; y++)        {            const unsigned char* src_prev = src_color;            const float* prev_factor = right_pass_factor;            const float* prev_color = right_pass_color;            // process 1st pixel separately since it has no previous            *right_pass_factor-- = 1.f;            for (int c = 0; c < channel; c++)            {                *right_pass_color-- = *src_color--;            }            // handle other pixels            for (int x = 1; x < width; x++)            {                // determine difference in pixel color between current and previous                // calculation is different depending on number of channels                int diff = getDiffFactor(src_color, src_color - 3);                //    src_prev = src_color;                float alpha_f = range_table_f[diff];                *right_pass_factor-- = inv_alpha_f + alpha_f * (*prev_factor--);                for (int c = 0; c < channel; c++)                {                    *right_pass_color-- = inv_alpha_f * (*src_color--) + alpha_f * (*prev_color--);                }            }        }    }    // vertical pass will be applied on top on horizontal pass, while using pixel differences from original image    // result color stored in 'm_left_pass_color' and vertical pass will use it as source color    {        float* img_out = m_left_pass_color; // use as temporary buffer        const float* left_pass_color = m_left_pass_color;        const float* left_pass_factor = m_left_pass_factor;        const float* right_pass_color = m_right_pass_color;        const float* right_pass_factor = m_right_pass_factor;        int width_height = width * height;        for (int i = 0; i < width_height; i++)        {            // average color divided by average factor            float factor = 1.f / ((*left_pass_factor++) + (*right_pass_factor++));            for (int c = 0; c < channel; c++)            {                *img_out++ = (factor * ((*left_pass_color++) + (*right_pass_color++)));            }        }    }    ///////////////    // Down pass    {        const float* src_color_hor = m_left_pass_color; // result of horizontal pass filter        const unsigned char* src_color = img_src;        float* down_pass_color = m_down_pass_color;        float* down_pass_factor = m_down_pass_factor;        const unsigned char* src_prev = src_color;        const float* prev_color = down_pass_color;        const float* prev_factor = down_pass_factor;        // 1st line done separately because no previous line        for (int x = 0; x < width; x++)        {            *down_pass_factor++ = 1.f;            for (int c = 0; c < channel; c++)            {                *down_pass_color++ = *src_color_hor++;            }            src_color += channel;        }        // handle other lines        for (int y = 1; y < height; y++)        {            for (int x = 0; x < width; x++)            {                // determine difference in pixel color between current and previous                // calculation is different depending on number of channels                int diff = getDiffFactor(src_color, src_prev);                src_prev += channel;                src_color += channel;                float alpha_f = range_table_f[diff];                *down_pass_factor++ = inv_alpha_f + alpha_f * (*prev_factor++);                for (int c = 0; c < channel; c++)                {                    *down_pass_color++ = inv_alpha_f * (*src_color_hor++) + alpha_f * (*prev_color++);                }            }        }    }    ///////////////    // Up pass    {        // start from end and then go up to begining         int last_index = width * height * channel - 1;        const unsigned char* src_color = img_src + last_index;        const float* src_color_hor = m_left_pass_color + last_index; // result of horizontal pass filter        float* up_pass_color = m_up_pass_color + last_index;        float* up_pass_factor = m_up_pass_factor + (width * height - 1);        //    const unsigned char* src_prev = src_color;        const float* prev_color = up_pass_color;        const float* prev_factor = up_pass_factor;        // 1st line done separately because no previous line        for (int x = 0; x < width; x++)        {            *up_pass_factor-- = 1.f;            for (int c = 0; c < channel; c++)            {                *up_pass_color-- = *src_color_hor--;            }            src_color -= channel;        }        // handle other lines        for (int y = 1; y < height; y++)        {            for (int x = 0; x < width; x++)            {                // determine difference in pixel color between current and previous                // calculation is different depending on number of channels                src_color -= channel;                int diff = getDiffFactor(src_color, src_color + width * channel);                float alpha_f = range_table_f[diff];                *up_pass_factor-- = inv_alpha_f + alpha_f * (*prev_factor--);                for (int c = 0; c < channel; c++)                {                    *up_pass_color-- = inv_alpha_f * (*src_color_hor--) + alpha_f * (*prev_color--);                }            }        }    }    ///////////////    // average result of vertical pass is written to output buffer    {        const float* down_pass_color = m_down_pass_color;        const float* down_pass_factor = m_down_pass_factor;        const float* up_pass_color = m_up_pass_color;        const float* up_pass_factor = m_up_pass_factor;        int width_height = width * height;        for (int i = 0; i < width_height; i++)        {            // average color divided by average factor            float factor = 1.f / ((*up_pass_factor++) + (*down_pass_factor++));            for (int c = 0; c < channel; c++)            {                *img_dst++ = (unsigned char)(factor * ((*up_pass_color++) + (*down_pass_color++)));            }        }    }}

这份代码把算法逻辑梳理得非常清晰。

一开始我也觉得这个代码效果和速度都很不错了。

但是没有仔细看,经过朋友的提醒,说是这个算法有严重的bug,会有线条瑕疵。

把对应的参数值调大之后,的确。特别明显的条纹效应,直接毁图了。

后来我经过一段时间的思考和研究,明确了算法的问题所在。

针对这个算法的修复,已经有明确思路了。

在这里,暂时先放出我经过稍加整理了OP_RBF的完整C代码版本。

做了一些代码逻辑上的优化,降低内存的使用。

 int   getDiffFactor(const unsigned char* color1, const unsigned char* color2, const  int & channels){    int final_diff;    int component_diff[4];    // find absolute difference between each component    for (int i = 0; i < channels; i++)    {        component_diff[i] = abs(color1[i] - color2[i]);    }    // based on number of components, produce a single difference value in the 0-255 range    switch (channels)    {    case 1:        final_diff = component_diff[0];        break;    case 2:        final_diff = ((component_diff[0] + component_diff[1]) >> 1);        break;    case 3:        final_diff = ((component_diff[0] + component_diff[2]) >> 2) + (component_diff[1] >> 1);        break;    case 4:        final_diff = ((component_diff[0] + component_diff[1] + component_diff[2] + component_diff[3]) >> 2);        break;    default:        final_diff = 0;    }    _ASSERT(final_diff >= 0 && final_diff <= 255);    return final_diff;}void CRB_HorizontalFilter(unsigned char* Input, unsigned char* Output, int Width, int Height, int Channels, float * range_table_f, float inv_alpha_f, float* left_Color_Buffer, float* left_Factor_Buffer, float* right_Color_Buffer, float* right_Factor_Buffer){    // Left pass and Right pass     int Stride = Width * Channels;    const unsigned char* src_left_color = Input;    float* left_Color = left_Color_Buffer;    float* left_Factor = left_Factor_Buffer;    int last_index = Stride * Height - 1;    const unsigned char* src_right_color = Input + last_index;    float* right_Color = right_Color_Buffer + last_index;    float* right_Factor = right_Factor_Buffer + Width * Height - 1;    for (int y = 0; y < Height; y++)    {        const unsigned char* src_left_prev = Input;        const float* left_prev_factor = left_Factor;        const float* left_prev_color = left_Color;        const unsigned char* src_right_prev = src_right_color;        const float* right_prev_factor = right_Factor;        const float* right_prev_color = right_Color;        // process 1st pixel separately since it has no previous        {            //if x = 0             *left_Factor++ = 1.f;            *right_Factor-- = 1.f;            for (int c = 0; c < Channels; c++)            {                *left_Color++ = *src_left_color++;                *right_Color-- = *src_right_color--;            }        }        // handle other pixels        for (int x = 1; x < Width; x++)        {            // determine difference in pixel color between current and previous            // calculation is different depending on number of channels            int left_diff = getDiffFactor(src_left_color, src_left_prev, Channels);            src_left_prev = src_left_color;            int right_diff = getDiffFactor(src_right_color, src_right_color - Channels, Channels);            src_right_prev = src_right_color;            float left_alpha_f = range_table_f[left_diff];            float right_alpha_f = range_table_f[right_diff];            *left_Factor++ = inv_alpha_f + left_alpha_f * (*left_prev_factor++);            *right_Factor-- = inv_alpha_f + right_alpha_f * (*right_prev_factor--);            for (int c = 0; c < Channels; c++)            {                *left_Color++ = (inv_alpha_f * (*src_left_color++) + left_alpha_f * (*left_prev_color++));                *right_Color-- = (inv_alpha_f * (*src_right_color--) + right_alpha_f * (*right_prev_color--));            }        }    }    // vertical pass will be applied on top on horizontal pass, while using pixel differences from original image    // result color stored in 'leftcolor' and vertical pass will use it as source color    {        unsigned char* dst_color = Output; // use as temporary buffer          const float* leftcolor = left_Color_Buffer;        const float* leftfactor = left_Factor_Buffer;        const float* rightcolor = right_Color_Buffer;        const float* rightfactor = right_Factor_Buffer;        int width_height = Width * Height;        for (int i = 0; i < width_height; i++)        {            // average color divided by average factor            float factor = 1.f / ((*leftfactor++) + (*rightfactor++));            for (int c = 0; c < Channels; c++)            {                *dst_color++ = (factor * ((*leftcolor++) + (*rightcolor++)));            }        }    }}void CRB_VerticalFilter(unsigned char* Input, unsigned char* Output, int Width, int Height, int Channels, float * range_table_f, float inv_alpha_f, float* down_Color_Buffer, float* down_Factor_Buffer, float* up_Color_Buffer, float* up_Factor_Buffer){    // Down pass and Up pass     int Stride = Width * Channels;    const unsigned char* src_color_first_hor = Output; // result of horizontal pass filter     const unsigned char* src_down_color = Input;    float* down_color = down_Color_Buffer;    float* down_factor = down_Factor_Buffer;    const unsigned char* src_down_prev = src_down_color;    const float* down_prev_color = down_color;    const float* down_prev_factor = down_factor;    int last_index = Stride * Height - 1;    const unsigned char* src_up_color = Input + last_index;    const unsigned char* src_color_last_hor = Output + last_index; // result of horizontal pass filter    float* up_color = up_Color_Buffer + last_index;    float* up_factor = up_Factor_Buffer + (Width * Height - 1);    const float* up_prev_color = up_color;    const float* up_prev_factor = up_factor;    // 1st line done separately because no previous line    {        //if y=0         for (int x = 0; x < Width; x++)        {            *down_factor++ = 1.f;            *up_factor-- = 1.f;            for (int c = 0; c < Channels; c++)            {                *down_color++ = *src_color_first_hor++;                *up_color-- = *src_color_last_hor--;            }            src_down_color += Channels;            src_up_color -= Channels;        }    }    // handle other lines     for (int y = 1; y < Height; y++)    {        for (int x = 0; x < Width; x++)        {            // determine difference in pixel color between current and previous            // calculation is different depending on number of channels            int down_diff = getDiffFactor(src_down_color, src_down_prev, Channels);            src_down_prev += Channels;            src_down_color += Channels;            src_up_color -= Channels;            int up_diff = getDiffFactor(src_up_color, src_up_color + Stride, Channels);            float down_alpha_f = range_table_f[down_diff];            float up_alpha_f = range_table_f[up_diff];            *down_factor++ = inv_alpha_f + down_alpha_f * (*down_prev_factor++);            *up_factor-- = inv_alpha_f + up_alpha_f * (*up_prev_factor--);            for (int c = 0; c < Channels; c++)            {                *down_color++ = inv_alpha_f * (*src_color_first_hor++) + down_alpha_f * (*down_prev_color++);                *up_color-- = inv_alpha_f * (*src_color_last_hor--) + up_alpha_f * (*up_prev_color--);            }        }    }    // average result of vertical pass is written to output buffer    {        unsigned char *dst_color = Output;        const float* downcolor = down_Color_Buffer;        const float* downfactor = down_Factor_Buffer;        const float* upcolor = up_Color_Buffer;        const float* upfactor = up_Factor_Buffer;        int width_height = Width * Height;        for (int i = 0; i < width_height; i++)        {            // average color divided by average factor            float factor = 1.f / ((*upfactor++) + (*downfactor++));            for (int c = 0; c < Channels; c++)            {                *dst_color++ = (factor * ((*upcolor++) + (*downcolor++)));            }        }    }}// memory must be reserved before calling image filter// this implementation of filter uses plain C++, single threaded// channel count must be 3 or 4 (alpha not used)void    CRBFilter(unsigned char* Input, unsigned char* Output, int Width, int Height, int Stride, float sigmaSpatial, float sigmaRange){    int Channels = Stride / Width;    int reserveWidth = Width;    int reserveHeight = Height;    // basic sanity check    _ASSERT(Input);    _ASSERT(Output);    _ASSERT(reserveWidth >= 10 && reserveWidth < 10000);    _ASSERT(reserveHeight >= 10 && reserveHeight < 10000);    _ASSERT(Channels >= 1 && Channels <= 4);    int reservePixels = reserveWidth * reserveHeight;    int numberOfPixels = reservePixels * Channels;    float* leftColorBuffer = (float*)calloc(sizeof(float)*numberOfPixels, 1);    float* leftFactorBuffer = (float*)calloc(sizeof(float)*reservePixels, 1);    float* rightColorBuffer = (float*)calloc(sizeof(float)*numberOfPixels, 1);    float* rightFactorBuffer = (float*)calloc(sizeof(float)*reservePixels, 1);    if (leftColorBuffer == NULL || leftFactorBuffer == NULL || rightColorBuffer == NULL || rightFactorBuffer == NULL)    {        if (leftColorBuffer)  free(leftColorBuffer);        if (leftFactorBuffer) free(leftFactorBuffer);        if (rightColorBuffer) free(rightColorBuffer);        if (rightFactorBuffer) free(rightFactorBuffer);        return;    }    float* downColorBuffer = leftColorBuffer;    float* downFactorBuffer = leftFactorBuffer;    float* upColorBuffer = rightColorBuffer;    float* upFactorBuffer = rightFactorBuffer;    // compute a lookup table    float alpha_f = static_cast<float>(exp(-sqrt(2.0) / (sigmaSpatial * 255)));    float inv_alpha_f = 1.f - alpha_f;    float range_table_f[255 + 1];    float inv_sigma_range = 1.0f / (sigmaRange * 255);    float ii = 0.f;    for (int i = 0; i <= 255; i++, ii -= 1.f)    {        range_table_f[i] = alpha_f * exp(ii * inv_sigma_range);    }    CRB_HorizontalFilter(Input, Output, Width, Height, Channels, range_table_f, inv_alpha_f, leftColorBuffer, leftFactorBuffer, rightColorBuffer, rightFactorBuffer);    CRB_VerticalFilter(Input, Output, Width, Height, Channels, range_table_f, inv_alpha_f, downColorBuffer, downFactorBuffer, upColorBuffer, upFactorBuffer);    if (leftColorBuffer)    {        free(leftColorBuffer);        leftColorBuffer = NULL;    }    if (leftFactorBuffer)    {        free(leftFactorBuffer);        leftFactorBuffer = NULL;    }    if (rightColorBuffer)    {        free(rightColorBuffer);        rightColorBuffer = NULL;    }    if (rightFactorBuffer)    {        free(rightFactorBuffer);        rightFactorBuffer = NULL;    }}

void    CRBFilter(unsigned char* Input, unsigned char* Output, int Width, int Height, int Stride, float sigmaSpatial, float sigmaRange);

就这样一个简单的接口函数。

代码也比较清晰,不多做解释了。

上个效果演示图。

前后对比,注意看人物的眼袋部分,有明显的条纹效应。

好久没上来发文了,最近有好多网友,询问关于《票据OCR前预处理 (附Demo)》一些算法相关的细节。

说实话,这个代码写得有点久了,逻辑记不太清楚了。

待后面有时间,把代码梳理一下,再写篇博客分享完整逻辑代码。

俺的联系方式如下:

邮箱: gaozhihan@vip.qq.com

联系我时请说明来意,不然一律忽略,谢谢。 

若此博文能帮到您,欢迎扫码小额赞助。

微信:  

 

支付宝: 

相关文章

    无相关信息

微信扫一扫

第七城市微信公众平台