磨皮对于现在的图像处理软件，可以说是一项重要的功能，在天天P图，可牛，ps，美图秀秀等软件中随处可见，有可能即使你非常熟悉图像处理的算法，然而却不懂磨皮怎么实现。其实磨皮就是所谓的保边缘滤波，也就是说在图像处理领域只要是滤波算法都可以实现磨皮，只是效果好坏的区别，然而现在对于大部分，都要求具有保细节的功能，这边先给大家介绍一种算法：双指数保边缘滤波，对应的 外围文献为：《Bi-Exponential Edge-Preserving Smoother》

下面是这篇文献最重要的部分，算法整个过程分为三个步骤：

（1）水平方向递归：包括对原始图像的水平方向进行向前递归（公式1）、对原始图片的水平方向进行向后递归（公式3），然后把向前递归结果、向后递归结果、原图像数据做一个加权组合(公式5)，得到新的像素值。

（2）垂直方向递归：同样的对原始图像数据进行与水平方向的递归方式类似的方法，计算得到新的像素值

（3）把垂直递归与水平递归的结果相加，然后除以2，得到最后的结果

这个算法具有高度并行的效果，因此启用多线程毫无压力，而且算法是图像中常用的递归加速，因此速度挺快的，如果还觉得速度不够快，可以采用查询表的方式

接着贴一下部分重要函数的代码，供参考，下面代码中我用了查询表的方式进行加速，速度可以比之前提高一倍，效果和原算法肉眼看不出来有什么区别，差别非常小，因此建议用查询表进行加速，还有其它的一些小细节我也没有根据原文进行写代码。

[cpp] view plaincopy
  1. void CBeeps::Run(BYTE* pImage, int nWidth, int nHeight, int nStride)
  2. {
  3. if (pImage == NULL)
  4. {
  5. return;
  6. }
  7. m_pImage=pImage;
  8. m_nWidth=nWidth;
  9. m_nHeight=nHeight;
  10. m_nStride=nStride;
  11. ApplyBiExponentialEdgePreservingSmoother(20,0.02);
  12. }
  13. void CBeeps::ApplyBiExponentialEdgePreservingSmoother(double photometricStandardDeviation, double spatialDecay)
  14. {
  15. m_exp_table=new double[256];
  16. m_g_table=new double[256];
  17. double c=-0.5/(photometricStandardDeviation * photometricStandardDeviation);
  18. double mu=spatialDecay/(2-spatialDecay);
  19. for (int i=0;i<=255;i++)
  20. {
  21. float a=exp(c*i*i);
  22. m_exp_table[i]=(1-spatialDecay)* exp(c*i*i);
  23. m_g_table[i]=mu*i;
  24. }
  25. BYTE* p0 =m_pImage;
  26. const int nChannel = 4;
  27. int m_length =m_nWidth*m_nHeight;
  28. float maxerror=0;
  29. float sum=0;
  30. // 对每个channel进行处理
  31. #pragma omp parallel for
  32. for (int idxChannel=0;idxChannel <nChannel; idxChannel++)
  33. {
  34. double *data1 = new double[m_length];
  35. double* data2 = new double[m_length];
  36. //double* data3 = new double[m_length];
  37. //double* data4 = new double[m_length];
  38. BYTE *p1=p0+idxChannel;
  39. for (int i = 0; i < m_length;++i)
  40. {
  41. data1[i] = p1[i * nChannel];
  42. }
  43. memcpy(data2,data1, sizeof(double) * m_length);
  44. //memcpy(data3,data1, sizeof(double) * m_length);
  45. //memcpy(data4,data1, sizeof(double) * m_length);
  46. //CBEEPSHorizontalVertical hv(data3, m_nWidth, m_nHeight,  photometricStandardDeviation, spatialDecay);
  47. //hv.run();
  48. //CBEEPSVerticalHorizontal vh(data4, m_nWidth, m_nHeight, photometricStandardDeviation, spatialDecay);
  49. //vh.run();
  50. runHorizontalVertical(data1, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
  51. runVerticalHorizontal(data2, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
  52. sum=0;
  53. for (int i =0;i<m_length;++i)
  54. {
  55. //double val = (data3[i] + data4[i]) * 0.5;
  56. double val=(data1[i] + data2[i]) * 0.5;
  57. //double error0=abs((int)val2-(int)val);
  58. //sum=sum+error0;
  59. //if (error0>maxerror)
  60. //{
  61. //maxerror=error0;
  62. //}
  63. if(255.0<val)val=255.0;
  64. p1[i * nChannel]=(BYTE)val;
  65. }
  66. //sum=sum/m_length;
  67. delete data1;
  68. delete data2;
  69. }//for
  70. delete m_exp_table;m_exp_table=NULL;
  71. delete m_g_table;m_g_table=NULL;
  72. }
  73. //垂直方向递归
  74. void CBeeps::runVerticalHorizontal(double *data,int width,int height,double spatialDecay,double *exp_table,double *g_table)
  75. {
  76. int length0=height*width;
  77. double* g= new double[length0];
  78. int m = 0;
  79. for (int k2 = 0;k2<height;++k2)
  80. {
  81. int n = k2;
  82. for (int k1 = 0;k1<width;++k1)
  83. {
  84. g[n]=data[m++];
  85. n += height;
  86. }
  87. }
  88. double*p = new double[length0];
  89. double*r = new double[length0];
  90. memcpy(p, g, sizeof(double) * length0);
  91. memcpy(r, g, sizeof(double) * length0);
  92. for (int k1 = 0;k1<width; ++k1)
  93. {
  94. int startIndex=k1 * height;
  95. double mu = 0.0;
  96. for (int k=startIndex+1,K =startIndex+height;k<K;++k)
  97. {
  98. int div0=fabs(p[k]-p[k-1]);
  99. mu =exp_table[div0];
  100. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文献中的公式1，这里做了一下修改，效果影响不大
  101. }
  102. for (int k =startIndex+height-2;startIndex <= k;--k)
  103. {
  104. int div0=fabs(r[k]-r[k+1]);
  105. mu =exp_table[div0];
  106. r[k] = r[k+1] * mu + r[k] * (1.0-mu) ;//文献公式3
  107. }
  108. }
  109. double rho0=1.0/(2-spatialDecay);
  110. for (int k = 0;k <length0;++k)
  111. {
  112. r[k]= (r[k]+p[k])*rho0-g_table[(int)g[k]];
  113. }
  114. m = 0;
  115. for (int k1=0;k1<width;++k1)
  116. {
  117. int n = k1;
  118. for (int k2 =0;k2<height;++k2)
  119. {
  120. data[n] = r[m++];
  121. n += width;
  122. }
  123. }
  124. memcpy(p,data, sizeof(double) * length0);
  125. memcpy(r,data, sizeof(double) * length0);
  126. for (int k2 = 0; k2<height;++k2)
  127. {
  128. int startIndex=k2 * width;
  129. double mu = 0.0;
  130. for (int k=startIndex+1, K=startIndex+width;k<K;++k)
  131. {
  132. int div0=fabs(p[k]-p[k-1]);
  133. mu =exp_table[div0];
  134. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
  135. }
  136. for (int k=startIndex+width-2;startIndex<=k;--k)
  137. {
  138. int div0=fabs(r[k]-r[k+1]);
  139. mu =exp_table[div0];
  140. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu) ;
  141. }
  142. }
  143. double init_gain_mu=spatialDecay/(2-spatialDecay);
  144. for (int k = 0; k <length0; k++)
  145. {
  146. data[k]=(p[k]+r[k])*rho0-data[k]*init_gain_mu;//文献中的公式5
  147. }
  148. delete p;
  149. delete r;
  150. delete g;
  151. }
  152. //水平方向递归
  153. void CBeeps::runHorizontalVertical(double *data,int width,int height,double spatialDecay,double *exptable,double *g_table)
  154. {
  155. int length=width*height;
  156. double* g = new double[length];
  157. double* p = new double[length];
  158. double* r = new double[length];
  159. memcpy(p,data, sizeof(double) * length);
  160. memcpy(r,data, sizeof(double) * length);
  161. double rho0=1.0/(2-spatialDecay);
  162. for (int k2 = 0;k2 < height;++k2)
  163. {
  164. int startIndex=k2 * width;
  165. for (int k=startIndex+1,K=startIndex+width;k<K;++k)
  166. {
  167. int div0=fabs(p[k]-p[k-1]);
  168. double mu =exptable[div0];
  169. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文献公式1
  170. }
  171. for (int k =startIndex + width - 2;startIndex <= k;--k)
  172. {
  173. int div0=fabs(r[k]-r[k+1]);
  174. double mu =exptable[div0];
  175. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);//文献公式3
  176. }
  177. for (int k =startIndex,K=startIndex+width;k<K;k++)
  178. {
  179. r[k]=(r[k]+p[k])*rho0- g_table[(int)data[k]];
  180. }
  181. }
  182. int m = 0;
  183. for (int k2=0;k2<height;k2++)
  184. {
  185. int n = k2;
  186. for (int k1=0;k1<width;k1++)
  187. {
  188. g[n] = r[m++];
  189. n += height;
  190. }
  191. }
  192. memcpy(p, g, sizeof(double) * height * width);
  193. memcpy(r, g, sizeof(double) * height * width);
  194. for (int k1=0;k1<width;++k1)
  195. {
  196. int startIndex=k1 * height;
  197. double mu = 0.0;
  198. for (int k =startIndex+1,K =startIndex+height;k<K;++k)
  199. {
  200. int div0=fabs(p[k]-p[k-1]);
  201. mu =exptable[div0];
  202. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
  203. }
  204. for (int k=startIndex+height-2;startIndex<=k;--k)
  205. {
  206. int div0=fabs(r[k]-r[k+1]);
  207. mu =exptable[div0];
  208. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);
  209. }
  210. }
  211. double init_gain_mu=spatialDecay/(2-spatialDecay);
  212. for (int k = 0;k <length;++k)
  213. {
  214. r[k]= (r[k]+p[k])*rho0- g[k]*init_gain_mu;
  215. }
  216. m = 0;
  217. for (int k1=0;k1<width;++k1)
  218. {
  219. int n = k1;
  220. for (int k2=0;k2<height;++k2)
  221. {
  222. data[n]=r[m++];
  223. n += width;
  224. }
  225. }
  226. delete p;
  227. delete r;
  228. delete g;
  229. }

最后看一下我用这个算法，实现磨皮的结果：

原图

美图秀秀磨皮结果

双指数递归算法磨皮结果

本文地址：http://blog.csdn.net/hjimce/article/details/45420965    作者：hjimce     联系qq：1393852684     更多资源请关注我的博客：http://blog.csdn.net/hjimce                  原创文章，版权所有，转载请保留这两行作者信息

图像处理（五）双指数磨皮相关推荐

  1. 双指数边缘平滑滤波器用于磨皮算法的尝试

    本篇博文来自博主Imageshop,打赏或想要查阅更多内容可以移步至Imageshop. 转载自:https://www.cnblogs.com/Imageshop/p/3293300.html   ...

  2. 双指数边缘平滑滤波器用于磨皮算法的尝试。

    说起为什么会看到这个东西,那还真的绕一圈.首先在写<Single Image Haze Removal Using Dark Channel Prior>一文中图像去雾算法的原理.实现.效 ...

  3. 双指数核脉冲信号分析

    双指数核脉冲信号常见形态与成因简析 前言 信号分析 表达式分析 信号形状与产生原因分析 小结 前言 双指数信号是在核信号处理领域很常见的一种型号,这是因为在单指数信号通过带电容电路后,都会变为双指数信 ...

  4. Holt-Winters双指数平滑的java实现

    当数据有趋势时,简单的指数平滑并不会很好,这很不方便.[1]在这种情况下,以"双指数平滑"或"二阶指数平滑"的名称设计了几种方法,这是指数滤波器的两次递归应用, ...

  5. 温度补偿 matlab,基于传感器温度补偿方法的双指数函数模型的温度补偿算法设计...

    描述 0 引言 目前基于传感器的温度补偿方法主要分为模拟硬件设计和数字信号处理两种方法.模拟硬件通常采用PTAT 和CTAT 等技术来设计读出电路.数字信号处理方法通常包括线性拟合.最小二乘多项式 拟 ...

  6. 《OpenCv视觉之眼》Python图像处理五 :Opencv图像去噪处理之均值滤波、方框滤波、中值滤波和高斯滤波

    本专栏主要介绍如果通过OpenCv-Python进行图像处理,通过原理理解OpenCv-Python的函数处理原型,在具体情况中,针对不同的图像进行不同等级的.不同方法的处理,以达到对图像进行去噪.锐 ...

  7. [Python图像处理] 五.图像融合、加法运算及图像类型转换

    该系列文章是讲解Python OpenCV图像处理知识,前期主要讲解图像入门.OpenCV基础用法,中期讲解图像处理的各种算法,包括图像锐化算子.图像增强技术.图像分割等,后期结合深度学习研究图像识别 ...

  8. 图像处理五：python读取图片的几种方式

    一.读取图片方式 PIL.opencv.scikit-image: (1)PIL和Pillow只提供最基础的数字图像处理,功能有限: (2)opencv实际上是一个c++库,只是提供了python接口 ...

  9. 《模拟电子技术基础》课程笔记(五)——双极型三极管

    目录 双极性三极管 三极管的结构 三极管的放大作用 三极管中载流子的运动过程 三极管的电流分配关系 三极管的特性曲线: 三极管的主要参数 PNP型三极管 双极性三极管 双极性三极管又称三极管.晶体管. ...

最新文章

  1. myeclipse java可视化_使用MyEclipse可视化开发Hibernate实例
  2. 胖爷的vim实用手册 - 基础篇（打开、关闭、移动、搜索）
  3. Java NIO学习系列六：Java中的IO模型
  4. 中国科协发布 2021 开源创新榜，阿里巴巴 2 大开源社区、5 大开源项目上榜
  5. WebService到底是什？
  6. Git Flow分支管理
  7. 最大规模线上新基建项目拉开大幕！第127届广交会今天正式开展
  8. python matlibplot_python matlibplot绘制3D图形
  9. 遍历集合的两种方式：迭代器和增强型for循环
  10. 虚拟串口 服务器,ZNetCManager
  11. 【李宏毅2020 ML/DL】P53-55 Conditional Generation by RNN Attention Pointer Network Recursive
  12. 混合线性规划matlab,matlab求解混合的非线性规划软件说明
  13. Goldendict 1.5.0 VS2015 Qt 5.7 源代码编译
  14. MP3/MP4原理电路图下载全搜集
  15. python中range函数是什么意思_python中range什么意思
  16. Python实现批量修改图片名称并存入新文件夹
  17. DCT域图像水印技术
  18. maya xgen基础头发
  19. CSS漂亮搜索框代码
  20. 【MXNet学习16】在MXNet中使用Dropout

热门文章

  1. php如何从左往右轮播,js实现从左向右滑动式轮播图效果
  2. 转载：谢谢原作者： 块设备驱动实战基础篇二 （继续完善170行过滤驱动代码至200行）
  3. 实战SSM_O2O商铺_36【商品】商品列表之Dao+Service+Controller层的实现
  4. MyBatis-25MyBatis缓存配置【集成Redis】
  5. Quartz-Job 详解
  6. 学习笔记Hive（九）—— 实例：航空客户价值分析数据预处理
  7. input 只输入数字并限制最大输入长度
  8. python入门之控制结构顺序与选择结构_Python 入门之控制结构 - 顺序与选择结构——第1关：顺序结构...
  9. spring注解@service(service)括号中的service有什么用？
  10. 2021-01-10 Halcon初学者知识 【10】形状匹配 【二】模板的形状匹配