高斯滤波器是图像处理中经常用到的滤波器,其滤波核函数为:
为简单起见,这里省略了归一化因子。
由的可分离特性:
得:
其中为输入图像,为输出图像,为滤波模板半径。根据准则,通常使。
由上式可见,我们可以将二维高斯滤波分解为两次一维高斯滤波。
对于二维高斯滤波,设图像大小,高斯模板大小,处理每个像素点需要次操作,则算法复杂度。若使用一维高斯核先对图像逐行滤波,再对中间结果逐列滤波,则处理每个像素需要次操作,算法复杂度,随着滤波模板尺寸的增大,算法优势越明显。
程序:
#include "stdafx.h"
#include<stdlib.h>
#include<math.h>
//边界处理
int Edge(int i,int x,int Max)
{
int k=i+x;
if(k<0)k=-i;
else if(k>=Max) k=Max-i-1;
else k=x;
return k;
}
//二维高斯处理灰度图像
extern "C" _declspec(dllexport) void GaussFilterGray(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize*nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilter=0.0;
double ImageData=0.0;
//生成二维高斯滤波核
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
int nDis_x = i-nCenter;
int nDis_y = j-nCenter;
pdKernel[j+i*nWindowSize]=exp(-(nDis_x*nDis_x+nDis_y*nDis_y)*scale2X);
dSum += pdKernel[i*nWindowSize+j];
}
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
pdKernel[i*nWindowSize+j] /= dSum;
}
}
//逐像素处理
for(int inx=0,i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth; j++,inx++)
{
dFilter=0;
//邻域内加权平均
for(int n=0,x=-nCenter; x<=nCenter; x++)
{
int i_x=Edge(i,x,nHeight);
for(int y=-nCenter; y<=nCenter; y++,n++)
{
int j_y=Edge(j,y,nWidth);
int index=(i+i_x)*Stride+j+j_y;//邻域内像素在内存中的下标
ImageData=buffer[index];
dFilter+=ImageData*pdKernel[n];
}
}
A[inx]= max(min(255,dFilter),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//一维高斯处理灰度图像
extern "C" _declspec(dllexport) void GaussFilterGray1D(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilter=0.0;
double ImageData=0.0;
//生成一维高斯核
for(int i=0; i<nWindowSize; i++)
{
int nDis_x = i-nCenter;
pdKernel[i]=exp(-(nDis_x*nDis_x)*scale2X);
dSum += pdKernel[i];
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
pdKernel[i] /= dSum;
}
//横向滤波
for(int inx=0,i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth; j++,inx++)
{
dFilter=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nWidth);
int index=inx+j_x;
ImageData=A[index];//从原图像A中取值
dFilter+=ImageData*pdKernel[n];
}
buffer[inx]= max(min(255,dFilter),0);//中间结果放在buffer中
}
}
//纵向滤波
for(int i=0;i<nWidth;i++)
{
for(int j=0;j<nHeight;j++)
{
dFilter=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nHeight);
int index=(j+j_x)*Stride+i;
ImageData=buffer[index];//从中间图像buffer中取值
dFilter+=ImageData*pdKernel[n];
}
A[j*Stride+i]= max(min(255,dFilter),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//二维高斯处理彩色图像
extern "C" _declspec(dllexport) void GaussFilterColor(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
int Step=3;
if(Stride==4*nWidth)Step=4;//四通道图像
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize*nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilterB=0.0;
double dFilterG=0.0;
double dFilterR=0.0;
int index;
double ImageData;
//生成二维高斯滤波核
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
int nDis_x = i-nCenter;
int nDis_y = j-nCenter;
pdKernel[j+i*nWindowSize]=exp(-(nDis_x*nDis_x+nDis_y*nDis_y)*scale2X);
dSum += pdKernel[i*nWindowSize+j];
}
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
pdKernel[i*nWindowSize+j] /= dSum;
}
}
for(int i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++)
{
int i_x=Edge(i,x,nHeight);
for(int y=-nCenter; y<=nCenter; y++,n++)
{
int j_y=Edge(j,y,nWidth);
index=(i+i_x)*Stride+(j+j_y)*Step;
//三通道BGR,四通道BGRA
ImageData=buffer[index];
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterR+=ImageData * pdKernel[n];
}
}
index=i*Stride+j*Step;
A[index]=max(min(dFilterB,255),0);
A[index+1]=max(min(dFilterG,255),0);
A[index+2]=max(min(dFilterR,255),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//一维高斯处理彩色图像
extern "C" _declspec(dllexport) void GaussFilterColor1D(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
int Step=3;
if(Stride==4*nWidth)Step=4;//四通道图像
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize];
double scale2X = 0.5/(dSigma*dSigma);
double dSum = 0.0;
double dFilterB=0;
double dFilterG=0;
double dFilterR=0;
double ImageData;
int index;
//一维高斯核
for(int i=0; i<nWindowSize; i++)
{
int nDis_x = i-nCenter;
pdKernel[i]=exp(-(nDis_x*nDis_x)*scale2X);
dSum += pdKernel[i];
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
pdKernel[i] /= dSum;
}
//横向滤波
for(int i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nWidth);
index=i*Stride+(j+j_x)*Step;
ImageData=A[index];//从原图像A中取值
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=A[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=A[index];
dFilterR+=ImageData * pdKernel[n];
}
index=i*Stride+j*Step;
buffer[index]=max(min(dFilterB,255),0);//中间结果放在buffer中
buffer[index+1]=max(min(dFilterG,255),0);
buffer[index+2]=max(min(dFilterR,255),0);
}
}
//纵向滤波
for(int i=0;i<nWidth;i++)
{
for(int j=0;j<nHeight;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nHeight);
int index=(j+j_x)*Stride+i*Step;
ImageData=buffer[index];//从中间图像buffer中取值
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterR+=ImageData * pdKernel[n];
}
index=j*Stride+i*Step;
A[index]=max(min(dFilterB,255),0);
A[index+1]=max(min(dFilterG,255),0);
A[index+2]=max(min(dFilterR,255),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
演示结果:
上面对一幅512*512的彩色图像,基本的高斯算法耗时1469ms,而快速高斯耗时439ms。选取的,此时的滤波模板大小为。
为简单起见,这里省略了归一化因子。
由的可分离特性:
得:
其中为输入图像,为输出图像,为滤波模板半径。根据准则,通常使。
由上式可见,我们可以将二维高斯滤波分解为两次一维高斯滤波。
对于二维高斯滤波,设图像大小,高斯模板大小,处理每个像素点需要次操作,则算法复杂度。若使用一维高斯核先对图像逐行滤波,再对中间结果逐列滤波,则处理每个像素需要次操作,算法复杂度,随着滤波模板尺寸的增大,算法优势越明显。
程序:
#include "stdafx.h"
#include<stdlib.h>
#include<math.h>
//边界处理
int Edge(int i,int x,int Max)
{
int k=i+x;
if(k<0)k=-i;
else if(k>=Max) k=Max-i-1;
else k=x;
return k;
}
//二维高斯处理灰度图像
extern "C" _declspec(dllexport) void GaussFilterGray(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize*nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilter=0.0;
double ImageData=0.0;
//生成二维高斯滤波核
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
int nDis_x = i-nCenter;
int nDis_y = j-nCenter;
pdKernel[j+i*nWindowSize]=exp(-(nDis_x*nDis_x+nDis_y*nDis_y)*scale2X);
dSum += pdKernel[i*nWindowSize+j];
}
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
pdKernel[i*nWindowSize+j] /= dSum;
}
}
//逐像素处理
for(int inx=0,i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth; j++,inx++)
{
dFilter=0;
//邻域内加权平均
for(int n=0,x=-nCenter; x<=nCenter; x++)
{
int i_x=Edge(i,x,nHeight);
for(int y=-nCenter; y<=nCenter; y++,n++)
{
int j_y=Edge(j,y,nWidth);
int index=(i+i_x)*Stride+j+j_y;//邻域内像素在内存中的下标
ImageData=buffer[index];
dFilter+=ImageData*pdKernel[n];
}
}
A[inx]= max(min(255,dFilter),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//一维高斯处理灰度图像
extern "C" _declspec(dllexport) void GaussFilterGray1D(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilter=0.0;
double ImageData=0.0;
//生成一维高斯核
for(int i=0; i<nWindowSize; i++)
{
int nDis_x = i-nCenter;
pdKernel[i]=exp(-(nDis_x*nDis_x)*scale2X);
dSum += pdKernel[i];
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
pdKernel[i] /= dSum;
}
//横向滤波
for(int inx=0,i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth; j++,inx++)
{
dFilter=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nWidth);
int index=inx+j_x;
ImageData=A[index];//从原图像A中取值
dFilter+=ImageData*pdKernel[n];
}
buffer[inx]= max(min(255,dFilter),0);//中间结果放在buffer中
}
}
//纵向滤波
for(int i=0;i<nWidth;i++)
{
for(int j=0;j<nHeight;j++)
{
dFilter=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nHeight);
int index=(j+j_x)*Stride+i;
ImageData=buffer[index];//从中间图像buffer中取值
dFilter+=ImageData*pdKernel[n];
}
A[j*Stride+i]= max(min(255,dFilter),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//二维高斯处理彩色图像
extern "C" _declspec(dllexport) void GaussFilterColor(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
int Step=3;
if(Stride==4*nWidth)Step=4;//四通道图像
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize*nWindowSize];
double dSum = 0.0;
double scale2X = 0.5/(dSigma*dSigma);
double dFilterB=0.0;
double dFilterG=0.0;
double dFilterR=0.0;
int index;
double ImageData;
//生成二维高斯滤波核
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
int nDis_x = i-nCenter;
int nDis_y = j-nCenter;
pdKernel[j+i*nWindowSize]=exp(-(nDis_x*nDis_x+nDis_y*nDis_y)*scale2X);
dSum += pdKernel[i*nWindowSize+j];
}
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
for(int j=0; j<nWindowSize; j++)
{
pdKernel[i*nWindowSize+j] /= dSum;
}
}
for(int i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++)
{
int i_x=Edge(i,x,nHeight);
for(int y=-nCenter; y<=nCenter; y++,n++)
{
int j_y=Edge(j,y,nWidth);
index=(i+i_x)*Stride+(j+j_y)*Step;
//三通道BGR,四通道BGRA
ImageData=buffer[index];
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterR+=ImageData * pdKernel[n];
}
}
index=i*Stride+j*Step;
A[index]=max(min(dFilterB,255),0);
A[index+1]=max(min(dFilterG,255),0);
A[index+2]=max(min(dFilterR,255),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
//一维高斯处理彩色图像
extern "C" _declspec(dllexport) void GaussFilterColor1D(unsigned char *A, int nWidth, int nHeight, int Stride,double dSigma)
{
int Step=3;
if(Stride==4*nWidth)Step=4;//四通道图像
unsigned char *buffer=(unsigned char*)malloc(Stride*nHeight);
memcpy(buffer,A,Stride*nHeight);
int nWindowSize = (int)(1+2*ceil(3*dSigma));
int nCenter = (nWindowSize)/2;
double* pdKernel = new double[nWindowSize];
double scale2X = 0.5/(dSigma*dSigma);
double dSum = 0.0;
double dFilterB=0;
double dFilterG=0;
double dFilterR=0;
double ImageData;
int index;
//一维高斯核
for(int i=0; i<nWindowSize; i++)
{
int nDis_x = i-nCenter;
pdKernel[i]=exp(-(nDis_x*nDis_x)*scale2X);
dSum += pdKernel[i];
}
//归一化
for(int i=0; i<nWindowSize; i++)
{
pdKernel[i] /= dSum;
}
//横向滤波
for(int i=0; i<nHeight; i++)
{
for(int j=0; j<nWidth;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nWidth);
index=i*Stride+(j+j_x)*Step;
ImageData=A[index];//从原图像A中取值
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=A[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=A[index];
dFilterR+=ImageData * pdKernel[n];
}
index=i*Stride+j*Step;
buffer[index]=max(min(dFilterB,255),0);//中间结果放在buffer中
buffer[index+1]=max(min(dFilterG,255),0);
buffer[index+2]=max(min(dFilterR,255),0);
}
}
//纵向滤波
for(int i=0;i<nWidth;i++)
{
for(int j=0;j<nHeight;j++)
{
dFilterB=0;
dFilterG=0;
dFilterR=0;
for(int n=0,x=-nCenter; x<=nCenter; x++,n++)
{
int j_x=Edge(j,x,nHeight);
int index=(j+j_x)*Stride+i*Step;
ImageData=buffer[index];//从中间图像buffer中取值
dFilterB+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterG+=ImageData * pdKernel[n];
index+=1;
ImageData=buffer[index];
dFilterR+=ImageData * pdKernel[n];
}
index=j*Stride+i*Step;
A[index]=max(min(dFilterB,255),0);
A[index+1]=max(min(dFilterG,255),0);
A[index+2]=max(min(dFilterR,255),0);
}
}
delete[]pdKernel;
delete[]buffer;
}
演示结果:
上面对一幅512*512的彩色图像,基本的高斯算法耗时1469ms,而快速高斯耗时439ms。选取的,此时的滤波模板大小为。