圖像保邊濾波器——雙邊濾波算法

雙邊濾波算法

图像保边滤波器——双边滤波算法
  1. #include "string.h"
  2. #include "stdio.h"
  3. #include "stdlib.h"
  4. #include "math.h"
  5. #include"SoftSkin.h"
  6. //垂直方向遞歸
  7. void runVerticalHorizontal(double *data,int width,int height,double spatialDecay,double *exp_table,double *g_table)
  8. {
  9. int length0=height*width;
  10. double* g= new double[length0];
  11. int m = 0;
  12. for (int k2 = 0;k2
  13. {
  14. int n = k2;
  15. for (int k1 = 0;k1
  16. {
  17. g[n]=data[m++];
  18. n += height;
  19. }
  20. }
  21. double*p = new double[length0];
  22. double*r = new double[length0];
  23. memcpy(p, g, sizeof(double) * length0);
  24. memcpy(r, g, sizeof(double) * length0);
  25. for (int k1 = 0;k1
  26. {
  27. int startIndex=k1 * height;
  28. double mu = 0.0;
  29. for (int k=startIndex+1,K =startIndex+height;k
  30. {
  31. int div0=fabs(p[k]-p[k-1]);
  32. mu =exp_table[div0];
  33. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文獻中的公式1,這裡做了一下修改,效果影響不大
  34. }
  35. for (int k =startIndex+height-2;startIndex <= k;--k)
  36. {
  37. int div0=fabs(r[k]-r[k+1]);
  38. mu =exp_table[div0];
  39. r[k] = r[k+1] * mu + r[k] * (1.0-mu) ;//文獻公式3
  40. }
  41. }
  42. double rho0=1.0/(2-spatialDecay);
  43. for (int k = 0;k
  44. {
  45. r[k]= (r[k]+p[k])*rho0-g_table[(int)g[k]];
  46. }
  47. m = 0;
  48. for (int k1=0;k1
  49. {
  50. int n = k1;
  51. for (int k2 =0;k2
  52. {
  53. data[n] = r[m++];
  54. n += width;
  55. }
  56. }
  57. memcpy(p,data, sizeof(double) * length0);
  58. memcpy(r,data, sizeof(double) * length0);
  59. for (int k2 = 0; k2
  60. {
  61. int startIndex=k2 * width;
  62. double mu = 0.0;
  63. for (int k=startIndex+1, K=startIndex+width;k
  64. {
  65. int div0=fabs(p[k]-p[k-1]);
  66. mu =exp_table[div0];
  67. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
  68. }
  69. for (int k=startIndex+width-2;startIndex<=k;--k)
  70. {
  71. int div0=fabs(r[k]-r[k+1]);
  72. mu =exp_table[div0];
  73. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu) ;
  74. }
  75. }
  76. double init_gain_mu=spatialDecay/(2-spatialDecay);
  77. for (int k = 0; k
  78. {
  79. data[k]=(p[k]+r[k])*rho0-data[k]*init_gain_mu;//文獻中的公式5
  80. }
  81. delete p;
  82. delete r;
  83. delete g;
  84. }
  85. //水平方向遞歸
  86. void runHorizontalVertical(double *data,int width,int height,double spatialDecay,double *exptable,double *g_table)
  87. {
  88. int length=width*height;
  89. double* g = new double[length];
  90. double* p = new double[length];
  91. double* r = new double[length];
  92. memcpy(p,data, sizeof(double) * length);
  93. memcpy(r,data, sizeof(double) * length);
  94. double rho0=1.0/(2-spatialDecay);
  95. for (int k2 = 0;k2 < height;++k2)
  96. {
  97. int startIndex=k2 * width;
  98. for (int k=startIndex+1,K=startIndex+width;k
  99. {
  100. int div0=fabs(p[k]-p[k-1]);
  101. double mu =exptable[div0];
  102. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文獻公式1
  103. }
  104. for (int k =startIndex + width - 2;startIndex <= k;--k)
  105. {
  106. int div0=fabs(r[k]-r[k+1]);
  107. double mu =exptable[div0];
  108. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);//文獻公式3
  109. }
  110. for (int k =startIndex,K=startIndex+width;k
  111. {
  112. r[k]=(r[k]+p[k])*rho0- g_table[(int)data[k]];
  113. }
  114. }
  115. int m = 0;
  116. for (int k2=0;k2
  117. {
  118. int n = k2;
  119. for (int k1=0;k1
  120. {
  121. g[n] = r[m++];
  122. n += height;
  123. }
  124. }
  125. memcpy(p, g, sizeof(double) * height * width);
  126. memcpy(r, g, sizeof(double) * height * width);
  127. for (int k1=0;k1
  128. {
  129. int startIndex=k1 * height;
  130. double mu = 0.0;
  131. for (int k =startIndex+1,K =startIndex+height;k
  132. {
  133. int div0=fabs(p[k]-p[k-1]);
  134. mu =exptable[div0];
  135. p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
  136. }
  137. for (int k=startIndex+height-2;startIndex<=k;--k)
  138. {
  139. int div0=fabs(r[k]-r[k+1]);
  140. mu =exptable[div0];
  141. r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);
  142. }
  143. }
  144. double init_gain_mu=spatialDecay/(2-spatialDecay);
  145. for (int k = 0;k
  146. {
  147. r[k]= (r[k]+p[k])*rho0- g[k]*init_gain_mu;
  148. }
  149. m = 0;
  150. for (int k1=0;k1
  151. {
  152. int n = k1;
  153. for (int k2=0;k2
  154. {
  155. data[n]=r[m++];
  156. n += width;
  157. }
  158. }
  159. delete p;
  160. delete r;
  161. delete g;
  162. }
  163. void ApplyBiExponentialEdgePreservingSmoother(double photometricStandardDeviation, double spatialDecay, unsigned char* m_pImage, int m_nWidth, int m_nHeight,int m_nStride)
  164. {
  165. double m_exp_table[256];
  166. double m_g_table[256];
  167. memset(m_exp_table, 0, 256);
  168. memset(m_g_table, 0, 256);
  169. double c=-0.5/(photometricStandardDeviation * photometricStandardDeviation);
  170. double mu=spatialDecay/(2-spatialDecay);
  171. for (int i=0;i<=255;i++)
  172. {
  173. float a=exp(c*i*i);
  174. m_exp_table[i]=(1-spatialDecay)* exp(c*i*i);
  175. m_g_table[i]=mu*i;
  176. }
  177. unsigned char* p0 =m_pImage;
  178. const int nChannel = 4;
  179. int m_length =m_nWidth*m_nHeight;
  180. float maxerror=0;
  181. float sum=0;
  182. // 對每個channel進行處理
  183. for (int idxChannel=0;idxChannel
  184. {
  185. double *data1 = new double[m_length];
  186. double* data2 = new double[m_length];
  187. unsigned char *p1=p0+idxChannel;
  188. for (int i = 0; i < m_length;++i)
  189. {
  190. data1[i] = p1[i * nChannel];
  191. }
  192. memcpy(data2,data1, sizeof(double) * m_length);
  193. runHorizontalVertical(data1, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
  194. runVerticalHorizontal(data2, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
  195. sum=0;
  196. for (int i =0;i
  197. {
  198. double val=(data1[i] + data2[i]) * 0.5;
  199. if(255.0
  200. p1[i * nChannel]=(unsigned char)val;
  201. }
  202. delete data1;
  203. delete data2;
  204. }
  205. }
  206. void f_Bilateralfilter(unsigned char* pImage, int nWidth, int nHeight, int nStride, double std)
  207. {
  208. if (pImage == NULL || std == 0)
  209. {
  210. return;
  211. }
  212. ApplyBiExponentialEdgePreservingSmoother(std,std * 0.001,pImage, nWidth, nHeight, nStride);
  213. }

效果如下:

图像保边滤波器——双边滤波算法

原圖

图像保边滤波器——双边滤波算法

雙邊濾波效果圖(Delta=20)

遞歸雙邊濾波是一種快速雙邊濾波算法,這裡本人是對四通道處理的,大家可以在YUV顏色空間對Y通道處理,UV不變,這樣速度會更快一些。


分享到:


相關文章: