高斯函数图形
高斯滤波器是一类根据高斯函数的形状来选择权值的线性平滑滤波器。高斯平滑滤波器对于抑制服从正态分布的噪声非常有效。一维零均值高斯函数为:
g(x)=exp( -x^2/(2 sigma^2)
其中,高斯分布参数Sigma决定了高斯函数的宽度。对于图像处理来说,常用二维零均值离散高斯函数作平滑滤波器。
高斯函数具有五个重要的性质,这些性质使得它在早期图像处理中特别有用.这些性质表明,高斯平滑滤波器无论在空间域还是在频率域都是十分有效的低通滤波器,且在实际图像处理中得到了工程人员的有效使用.高斯函数具有五个十分重要的性质,它们是:
(1)二维高斯函数具有旋转对称性,即滤波器在各个方向上的平滑程度是相同的.一般来说,一幅图像的边缘方向是事先不知道的,因此,在滤波前是无法确定一个方向上比另一方向上需要更多的平滑.旋转对称性意味着高斯平滑滤波器在后续边缘检测中不会偏向任一方向.
(2)高斯函数是单值函数.这表明,高斯滤波器用像素邻域的加权均值来代替该点的像素值,而每一邻域像素点权值是随该点与中心点的距离单调增减的.这一性质是很重要的,因为边缘是一种图像局部特征,如果平滑运算对离算子中心很远的像素点仍然有很大作用,则平滑运算会使图像失真.
(3)高斯函数的傅立叶变换频谱是单瓣的.正如下面所示,这一性质是高斯函数付立叶变换等于高斯函数本身这一事实的直接推论.图像常被不希望的高频信号所污染(噪声和细纹理).而所希望的图像特征(如边缘),既含有低频分量,又含有高频分量.高斯函数付立叶变换的单瓣意味着平滑图像不会被不需要的高频信号所污染,同时保留了大部分所需信号.
(4)高斯滤波器宽度(决定着平滑程度)是由参数σ表征的,而且σ和平滑程度的关系是非常简单的.σ越大,高斯滤波器的频带就越宽,平滑程度就越好.通过调节平滑程度参数σ,可在图像特征过分模糊(过平滑)与平滑图像中由于噪声和细纹理所引起的过多的不希望突变量(欠平滑)之间取得折衷.
(5)由于高斯函数的可分离性,较大尺寸的高斯滤波器可以得以有效地实现.二维高斯函数卷积可以分两步来进行,首先将图像与一维高斯函数进行卷积,然后将卷积结果与方向垂直的相同一维高斯函数卷积.因此,二维高斯滤波的计算量随滤波模板宽度成线性增长而不是成平方增长.
g(x,y)={f(x-1,y-1)+f(x-1,y+1)+f(x+1,y-1)+f(x+1,y+1)+[f(x-1,y)+f(x,y-1)+f(x+1,y)+f(x,y+1)]*2+f(x,y)*4}/16;
其中,f(x,y)为图像中(x,y)点的灰度值,g(x,y)为该点经过高斯滤波后的值。
分三个部分:
行向量,列向量,矩阵类型
具体如下:
一维高斯kernel行向量:
/
//kernel X
Mat kernel_x;
int kernel_size = 5;
int radius = kernel_size/2;
double weight = 1.0;
kernel_x = Mat_<double>(1,kernel_size);
for (int i = -radius; i <= radius; i++)
kernel_x.at<double>(0,radius + i) = (double)exp(-((i * i)/(2 * weight * weight)) );
cout<<"kernel X vector = "<<kernel_x<<endl;
结果显示:
一维高斯kernel列向量:
/
//kernel Y
Mat kernel_y;
kernel_y = Mat_<double>(kernel_size,1);
for (int i = -radius; i <= radius; i++)
kernel_y.at<double>(radius + i,0) = (double)exp(-((i * i)/(2 * weight * weight)) );
cout<<"kernel Y vector = "<<endl<<kernel_y<<endl;
结果如下:
二维高斯kernel:
代码如下:
//kernel XY
Mat kernel_xy;
kernel_xy = Mat_<double>(kernel_size, kernel_size);
for (int i = -radius; i <= radius; i++){
for (int j = -radius; j <= radius; j++){
kernel_xy.at<double>(radius + i,radius + j) = exp((-((i * i) + (j * j))/(2 * weight * weight)));
}
}
cout<<"kernel XY vector = "<<endl<<kernel_xy<<endl;
结果如下:
完整高斯核模板产生代码如下:
#include <iostream>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;
int main()
{
/
//kernel X
Mat kernel_x;
int kernel_size = 5;
int radius = kernel_size/2;
double weight = 1.0;
kernel_x = Mat_<double>(1,kernel_size);
for (int i = -radius; i <= radius; i++)
kernel_x.at<double>(0,radius + i) = (double)exp(-((i * i)/(2 * weight * weight)) );
cout<<"kernel X vector = "<<kernel_x<<endl;
/
//kernel Y
Mat kernel_y;
kernel_y = Mat_<double>(kernel_size,1);
for (int i = -radius; i <= radius; i++)
kernel_y.at<double>(radius + i,0) = (double)exp(-((i * i)/(2 * weight * weight)) );
cout<<"kernel Y vector = "<<endl<<kernel_y<<endl;
//kernel XY
Mat kernel_xy;
kernel_xy = Mat_<double>(kernel_size, kernel_size);
for (int i = -radius; i <= radius; i++){
for (int j = -radius; j <= radius; j++){
kernel_xy.at<double>(radius + i,radius + j) = exp((-((i * i) + (j * j))/(2 * weight * weight)));
}
}
cout<<"kernel XY vector = "<<endl<<kernel_xy<<endl;
system("pause");
return 0;
}
——————————————————————————————————————
滤波实现可以如下:
#include "stdafx.h"
#include "highgui.h"
#include "cv.h"
//高斯滤波函数,在opencv里可以使用cvSmooth
void gaussianFilter(uchar* data, int width, int height)
{ int i, j, index, sum;
int templates[9] = { 1, 2, 1,
2, 4, 2,
1, 2, 1 };//模板的值
sum = height * width * sizeof(uchar);//图像所占内存的大小
uchar *tmpdata = (uchar*)malloc(sum);
memcpy((char*)tmpdata,(char*)data, sum);
for(i = 1;i < height - 1;i++)
{
for(j = 1;j < width - 1;j++)
{
index = sum = 0;
for(int m = i - 1;m < i + 2;m++)
{
for(int n = j - 1; n < j + 2;n++)
{
sum +=
tmpdata[m * width + n] *
templates[index++]; //处理
}
}
data[i * width + j] = sum / 16;
}
}
free(tmpdata);
}
void imgOperate( IplImage* image )
{
cvNamedWindow( "image-in", CV_WINDOW_AUTOSIZE );
cvNamedWindow( "image-out", CV_WINDOW_AUTOSIZE);
cvShowImage( "image-in", image );
/*
IplImage* out = cvCreateImage(
cvGetSize(image),
IPL_DEPTH_8U,
3
);
*/
//将色彩图像强制转化为灰度图像
IplImage* pGrayImg=NULL;
pGrayImg=cvCreateImage(cvGetSize(image),8,1);
cvCvtColor(image,pGrayImg,CV_RGB2GRAY);
// cvSmooth( image, out, CV_GAUSSIAN, 5,5 );
gaussianFilter((unsigned char*)pGrayImg->imageData,pGrayImg->width,pGrayImg->height);
cvShowImage( "image-out", pGrayImg );
cvReleaseImage( &pGrayImg );
cvWaitKey( 0 );
cvDestroyWindow("image-in" );
cvDestroyWindow("image-out" );
}
int main( int argc, char** argv )
{
IplImage* img = cvLoadImage("d:\\1.jpg");
imgOperate( img );
cvReleaseImage( &img );
return 0;
}
或者如下:
#include <iostream>
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
using namespace std;
using namespace cv;
void gaussianFilter(const Mat& src, Mat& result, int besarKernel, double delta);
void gaussianKernelGenerator(Mat &result, int besarKernel, double delta);
void rgb2gray(const Mat src, Mat &result);
int main(int argc, char *argv[])
{
Mat src = imread("lena.jpg");
rgb2gray(src, src);
Mat gausianFilter;
gaussianFilter(src, gausianFilter, 7, 1);
namedWindow("asli");
imshow("asli", src);
namedWindow("gaussian");
imshow("gaussian", gausianFilter);
cout<<gausianFilter.size()<<endl;
waitKey(0);
}
void gaussianFilter(const Mat& src, Mat& result, int besarKernel, double delta)
{
Mat kernel;
//inisialisasi kernel gaussians
gaussianKernelGenerator(kernel, besarKernel, delta);
int filterOffset = besarKernel / 2;
result = Mat::zeros(src.rows - filterOffset*2, src.cols - filterOffset*2, src.type());
double sumKeabuan;
double sumkernel;
for(int ysrc = filterOffset; ysrc < src.rows - filterOffset; ++ysrc){
for(int xsrc = filterOffset; xsrc < src.cols - filterOffset; ++xsrc){
sumKeabuan = 0;
sumkernel = 0;
for(int ykernel = -filterOffset; ykernel <= filterOffset; ++ykernel){
for(int xkernel = -filterOffset; xkernel <= filterOffset; ++xkernel){
sumKeabuan += src.at<uchar>(ysrc + ykernel, xsrc + xkernel) * kernel.at<double>(filterOffset + ykernel, filterOffset + xkernel);
sumkernel += kernel.at<double>(filterOffset + ykernel, filterOffset + xkernel);
}
}
result.at<uchar>(ysrc - filterOffset, xsrc - filterOffset) = sumKeabuan/sumkernel;
}
}
}
void gaussianKernelGenerator(Mat &result, int besarKernel, double delta)
{
int kernelRadius = besarKernel / 2;
result = Mat_<double>(besarKernel, besarKernel);
double pengali = 1 / (2 * (22 / 7) * delta * delta);
for(int filterY = - kernelRadius; filterY <= kernelRadius; filterY++){
for(int filterX = - kernelRadius; filterX <= kernelRadius; filterX++){
result.at<double>
(filterY + kernelRadius, filterX + kernelRadius) = exp(-( ( ( filterX * filterX ) + ( filterY * filterY ) ) / (delta * delta * 2) )) * pengali;
}
}
}
void rgb2gray(const Mat src, Mat &result)
{
CV_Assert(src.depth() != sizeof(uchar)); //harus 8 bit
result = Mat::zeros(src.rows, src.cols, CV_8UC1); //buat matrik 1 chanel
uchar data;
if(src.channels() == 3){
for( int i = 0; i < src.rows; ++i)
for( int j = 0; j < src.cols; ++j )
{
data = (uchar)(((Mat_<Vec3b>) src)(i,j)[0] * 0.0722 + ((Mat_<Vec3b>) src)(i,j)[1] * 0.7152 + ((Mat_<Vec3b>) src)(i,j)[2] * 0.2126);
result.at<uchar>(i,j) = data;
}
}else{
result = src;
}
}
但是方法二没有考虑边界的问题,结果图像比原图像小。
方法三:
方法三考虑了边界问题,并用reflect子函数进行了判断!
如题如下:
#include <iostream>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;
int reflect(int M, int x)
{
if(x < 0)
return (-x-1);
if(x >= M)
return (2*M-x-1);
else
return x;
}
double gaussianFunc(double x, double mu, double sigma)
{
return ((1/(sigma * CV_PI))*exp(-(((x-mu)/(sigma))*((x-mu)/(sigma)))/2.0)); //normal distribution function
}
vector<double> gaussianKernel(int kernelRadius, double sigma)
{
const int kernelSize = 2*kernelRadius+1;
vector<double> kernel(kernelSize); //kernel matrix
double sum = 0.0;
for(int i = -kernelRadius; i <= kernelRadius; i++)
{
kernel[i+kernelRadius] = gaussianFunc(i, kernelRadius, sigma);
sum += kernel[i+kernelRadius];
}
for(int i = 0; i < kernelSize; i++) //normalization
kernel[i] /= sum;
cout<<"kernel size = "<<kernel.size()<<endl;
return kernel;
}
void gaussianBlur(cv::Mat source, cv::Mat &target, double sigma, int kernelRadius)
{
const int kernelSize = 2*kernelRadius+1;
vector<double> kernel = gaussianKernel(kernelRadius, sigma);
float sum = 0.0;
float x1, y1;
target = source.clone();
cv::Mat temp = source.clone();
//along y-direction
for(int imgRow = 0; imgRow < source.rows; imgRow++)
{
for(int imgCol = 0; imgCol < source.cols; imgCol++)
{
sum = 0.0;
for(int i = -kernelRadius; i <= kernelRadius; i++)
{
y1 = reflect(source.rows, imgRow-i);
sum = sum + kernel[i+kernelRadius]*source.at<uchar>(y1, imgCol);
}
temp.at<uchar>(imgRow, imgCol) = sum;
}
}
//along x-direction
for(int imgRow = 0; imgRow < source.rows; imgRow++)
{
for(int imgCol = 0; imgCol < source.cols; imgCol++)
{
sum = 0.0;
for(int i = -kernelRadius; i <= kernelRadius; i++)
{
x1 = reflect(source.cols, imgCol-i);
sum += kernel[i+kernelRadius]*temp.at<uchar>(imgRow, x1);
}
target.at<uchar>(imgRow, imgCol) = sum;
}
}
}
int main()
{
Mat src = imread("lena.jpg",CV_LOAD_IMAGE_GRAYSCALE);
Mat dst;
imshow("src",src);
gaussianBlur(src, dst, 2.0, 9);
imshow("dst",dst);
cout<<dst.size()<<endl;
waitKey(0);
system("pause");
return 0;
}
结果显示如下:
如下方式四:
也考虑了边界问题,并且得到滤波结果:
代码如下:
#include <iostream>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;
int main()
{
cv::Mat image1 = cv::imread("lena.jpg",CV_LOAD_IMAGE_COLOR);
int kSize = 7;
float sigmax = 2;
float sigmay = 2;
int mean = kSize/2;
cv::Mat grey;
cv::Mat ggrey(image1.rows, image1.cols, CV_64F, cv::Scalar(100));
cv::Mat gauss(kSize, kSize, CV_64F, cv::Scalar(100));
cv::cvtColor(image1,grey,CV_BGR2GRAY);
grey.convertTo(grey,CV_64F);
double sum = 0.0;
for(int u = 0; u<kSize; u++)
{
for(int v = 0; v<kSize;v++)
{
gauss.at<double>(u,v) = (double)(exp(-0.5*((u-mean)*(u-mean)/(sigmax*sigmax)+(v-mean)*(v-mean)/(sigmay*sigmay)))/(2*3.141*sigmax*sigmay));
sum += gauss.at<double>(u,v);
}
}
for(int u = 0; u<kSize; ++u)
{
for(int v = 0; v<kSize;++v)
{
gauss.at<double>(u,v) = gauss.at<double>(u,v)/sum;
}
}
cout<<"gaussian = "<<gauss<<endl;
double d = 0.0, tmp = 0.0;
for(int i = 0; i<grey.rows; ++i)
{
for(int j = 0; j<grey.cols; ++j)
{
d = 0.0;
for(int u = -(kSize -1)/2; u<=(kSize - 1)/2;++u)
{
for(int v = -(kSize -1)/2; v<=(kSize - 1)/2;++v)
{
if((u+i)>=0 && (v+j)>=0 && (u+i)<grey.rows-1 && (v+j)<grey.cols-1)
{
d +=(grey.at<double>(i+u,j+v)/256.0)*gauss.at<double>(u+(kSize - 1)/2, v+(kSize - 1)/2);
}
}
}
ggrey.at<double>(i,j) = (double)(d); /// tmp);
}
}
/**/
grey.convertTo(grey,CV_8UC1);
cv::Mat gaussByCvLibrary(grey.size(),grey.type(),Scalar(0));
cv::GaussianBlur(grey,gaussByCvLibrary,cv::Size(kSize,kSize),sigmax);
cv::imshow("Library Gaussian", gaussByCvLibrary);
cv::imshow("My Gaussian", ggrey);
cv::waitKey(0);
system("pause");
return 0;
}
结果显示如下: