雙邊濾波算法
- #include "string.h"
- #include "stdio.h"
- #include "stdlib.h"
- #include "math.h"
- #include"SoftSkin.h"
- //垂直方向遞歸
- void runVerticalHorizontal(double *data,int width,int height,double spatialDecay,double *exp_table,double *g_table)
- {
- int length0=height*width;
- double* g= new double[length0];
- int m = 0;
- for (int k2 = 0;k2
- {
- int n = k2;
- for (int k1 = 0;k1
- {
- g[n]=data[m++];
- n += height;
- }
- }
- double*p = new double[length0];
- double*r = new double[length0];
- memcpy(p, g, sizeof(double) * length0);
- memcpy(r, g, sizeof(double) * length0);
- for (int k1 = 0;k1
- {
- int startIndex=k1 * height;
- double mu = 0.0;
- for (int k=startIndex+1,K =startIndex+height;k
- {
- int div0=fabs(p[k]-p[k-1]);
- mu =exp_table[div0];
- p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文獻中的公式1,這裡做了一下修改,效果影響不大
- }
- for (int k =startIndex+height-2;startIndex <= k;--k)
- {
- int div0=fabs(r[k]-r[k+1]);
- mu =exp_table[div0];
- r[k] = r[k+1] * mu + r[k] * (1.0-mu) ;//文獻公式3
- }
- }
- double rho0=1.0/(2-spatialDecay);
- for (int k = 0;k
- {
- r[k]= (r[k]+p[k])*rho0-g_table[(int)g[k]];
- }
- m = 0;
- for (int k1=0;k1
- {
- int n = k1;
- for (int k2 =0;k2
- {
- data[n] = r[m++];
- n += width;
- }
- }
- memcpy(p,data, sizeof(double) * length0);
- memcpy(r,data, sizeof(double) * length0);
- for (int k2 = 0; k2
- {
- int startIndex=k2 * width;
- double mu = 0.0;
- for (int k=startIndex+1, K=startIndex+width;k
- {
- int div0=fabs(p[k]-p[k-1]);
- mu =exp_table[div0];
- p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
- }
- for (int k=startIndex+width-2;startIndex<=k;--k)
- {
- int div0=fabs(r[k]-r[k+1]);
- mu =exp_table[div0];
- r[k] = r[k + 1] * mu + r[k] * (1.0 - mu) ;
- }
- }
- double init_gain_mu=spatialDecay/(2-spatialDecay);
- for (int k = 0; k
- {
- data[k]=(p[k]+r[k])*rho0-data[k]*init_gain_mu;//文獻中的公式5
- }
- delete p;
- delete r;
- delete g;
- }
- //水平方向遞歸
- void runHorizontalVertical(double *data,int width,int height,double spatialDecay,double *exptable,double *g_table)
- {
- int length=width*height;
- double* g = new double[length];
- double* p = new double[length];
- double* r = new double[length];
- memcpy(p,data, sizeof(double) * length);
- memcpy(r,data, sizeof(double) * length);
- double rho0=1.0/(2-spatialDecay);
- for (int k2 = 0;k2 < height;++k2)
- {
- int startIndex=k2 * width;
- for (int k=startIndex+1,K=startIndex+width;k
- {
- int div0=fabs(p[k]-p[k-1]);
- double mu =exptable[div0];
- p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);//文獻公式1
- }
- for (int k =startIndex + width - 2;startIndex <= k;--k)
- {
- int div0=fabs(r[k]-r[k+1]);
- double mu =exptable[div0];
- r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);//文獻公式3
- }
- for (int k =startIndex,K=startIndex+width;k
- {
- r[k]=(r[k]+p[k])*rho0- g_table[(int)data[k]];
- }
- }
- int m = 0;
- for (int k2=0;k2
- {
- int n = k2;
- for (int k1=0;k1
- {
- g[n] = r[m++];
- n += height;
- }
- }
- memcpy(p, g, sizeof(double) * height * width);
- memcpy(r, g, sizeof(double) * height * width);
- for (int k1=0;k1
- {
- int startIndex=k1 * height;
- double mu = 0.0;
- for (int k =startIndex+1,K =startIndex+height;k
- {
- int div0=fabs(p[k]-p[k-1]);
- mu =exptable[div0];
- p[k] = p[k - 1] * mu + p[k] * (1.0 - mu);
- }
- for (int k=startIndex+height-2;startIndex<=k;--k)
- {
- int div0=fabs(r[k]-r[k+1]);
- mu =exptable[div0];
- r[k] = r[k + 1] * mu + r[k] * (1.0 - mu);
- }
- }
- double init_gain_mu=spatialDecay/(2-spatialDecay);
- for (int k = 0;k
- {
- r[k]= (r[k]+p[k])*rho0- g[k]*init_gain_mu;
- }
- m = 0;
- for (int k1=0;k1
- {
- int n = k1;
- for (int k2=0;k2
- {
- data[n]=r[m++];
- n += width;
- }
- }
- delete p;
- delete r;
- delete g;
- }
- void ApplyBiExponentialEdgePreservingSmoother(double photometricStandardDeviation, double spatialDecay, unsigned char* m_pImage, int m_nWidth, int m_nHeight,int m_nStride)
- {
- double m_exp_table[256];
- double m_g_table[256];
- memset(m_exp_table, 0, 256);
- memset(m_g_table, 0, 256);
- double c=-0.5/(photometricStandardDeviation * photometricStandardDeviation);
- double mu=spatialDecay/(2-spatialDecay);
- for (int i=0;i<=255;i++)
- {
- float a=exp(c*i*i);
- m_exp_table[i]=(1-spatialDecay)* exp(c*i*i);
- m_g_table[i]=mu*i;
- }
- unsigned char* p0 =m_pImage;
- const int nChannel = 4;
- int m_length =m_nWidth*m_nHeight;
- float maxerror=0;
- float sum=0;
- // 對每個channel進行處理
- for (int idxChannel=0;idxChannel
- {
- double *data1 = new double[m_length];
- double* data2 = new double[m_length];
- unsigned char *p1=p0+idxChannel;
- for (int i = 0; i < m_length;++i)
- {
- data1[i] = p1[i * nChannel];
- }
- memcpy(data2,data1, sizeof(double) * m_length);
- runHorizontalVertical(data1, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
- runVerticalHorizontal(data2, m_nWidth, m_nHeight,spatialDecay,m_exp_table,m_g_table);
- sum=0;
- for (int i =0;i
- {
- double val=(data1[i] + data2[i]) * 0.5;
- if(255.0
- p1[i * nChannel]=(unsigned char)val;
- }
- delete data1;
- delete data2;
- }
- }
- void f_Bilateralfilter(unsigned char* pImage, int nWidth, int nHeight, int nStride, double std)
- {
- if (pImage == NULL || std == 0)
- {
- return;
- }
- ApplyBiExponentialEdgePreservingSmoother(std,std * 0.001,pImage, nWidth, nHeight, nStride);
- }
效果如下:
原圖
雙邊濾波效果圖(Delta=20)
遞歸雙邊濾波是一種快速雙邊濾波算法,這裡本人是對四通道處理的,大家可以在YUV顏色空間對Y通道處理,UV不變,這樣速度會更快一些。
閱讀更多 中國雲計算 的文章