高斯滤波,本文主要讲其如何通过
C
语言实现。不太擅长写理论性质的文章,这里仅仅阐述自己怎么实现以及简单的优化过程。
通常我们对获取的图像进行进一步处理时,往往需要先进行一次降噪,而通常我们选择的是高斯滤波也称高斯模糊,其因为较好的降噪效果而被广泛使用,具体的原因大家可以做实验以及看相关的论文去理解。
如果你经常使用MATLAB或者Python,可能只需要调用一下库函数即可实现。以matlab为例,其高斯滤波的常见使用,如下:
% Generate Kernel[5, 5], sigma = 1.0.
gaussKernel = fspecial('gaussian', [5, 5], 1.0);
% Compute gaussian filter result.
gaussResult = imfilter(Image, gaussKernel, 'replicate');
真的是方便至极,但自己不实现一下,好像没什么意思。并且嵌入式设备没有相应的库函数让你使用,所有还是有必要对其进行理解和学习。以下为正文分割线。
理论公式与模板生成
高斯滤波作为一种线性滤波器,其本质是带权值的加权均值滤波,相对于均值滤波,高斯滤波除了中心位置为1,其他位置根据距离的远近,越远滤波系数权重越低。
对高斯滤波的公式进行理解学习,参考如下公式:
其中,
为模板中心坐标,即坐标系原点。有了该公式,我们就可以实现高斯模板的生成了,参考如下
C
代码:
void kernelGenerate(float *gaussKernel, int k_size, float k_sigma) {
float flag = 1 / (2 * PI * k_sigma * k_sigma);
float coef;
float sum = 0;
int center = (k_size - 1) / 2;
for (int i = 0; i < k_size; i++) {
float x2 = (float)pow((i - center), 2);
for (int j = 0; j < k_size; j++) {
float y2 = (float)pow((j - center), 2);
coef = (-x2 - y2) * flag;
gaussKernel[i * k_size + j] = flag * exp(coef);
sum += gaussKernel[i * k_size + j];
}
}
/* Normalize */
if (sum) {
for (int i = 0; i < k_size * k_size; i++) {
gaussKernel[i] /= sum;
}
}
}
这个模板是生成的小数模板,要是生成整数模板呢?如下,将左上角置1。然后再进行取整操作。
void kernelGenerate(float *gaussKernel, int k_size, float k_sigma) {
float flag = 1 / (2 * PI * k_sigma * k_sigma);
float coef;
int center = (k_size - 1) / 2;
for (int i = 0; i < k_size; i++) {
float x2 = (float)pow((i - center), 2);
for (int j = 0; j < k_size; j++) {
float y2 = (float)pow((j - center), 2);
coef = (-x2 - y2) * flag;
gaussKernel[i * k_size + j] = flag * exp(coef);
}
}
/* Intergel Normalize */
if (gaussKernel[0]) {
float normal = 1.0 / gaussKernel[0];
for (int i = 0; i < k_size * k_size; i++) {
gaussKernel[i] *= normal;
}
}
}
有了模板之后,就需要进行滤波了。通过公式,可以发现滤波最关键的就是滤波的核的标准差,因此需要根据实际的需要对其进行相应的设置。
滤波处理
如果直接使用公式,进行滤波则可以写出如下代码:
void gaussian_filter(unsigned char *src, unsigned char *dst, int src_w, int src_h, int k_size, float k_sigma, int padd_type, int out_type) {
int img_w, img_h;
if (k_size == 1)
dst = src;
if (k_size % 2 == 0) {
printf("The kernel size not odd and ksize > 1! Convert it to k_size + 1");
k_size += 1;
}
int padd_size = (k_size - 1) / 2;
img_w = src_w + padd_size * 2;
img_h = src_h + padd_size * 2;
unsigned char *pBlur = (unsigned char *)malloc(img_w * img_h);
unsigned char *pDst = dst;
float *gaussKernel = (float *)malloc(sizeof(float) * k_size * k_size);
for (int i = 0; i < src_h; i++) {
unsigned char *pTmp = pBlur + (padd_size + i)* img_w + padd_size;
memcpy(pTmp, src + i * src_w, src_w);
/* col padding */
for (int j = 0; j < padd_size; j++) {
if (padd_type == 0) { // zero padding.
pTmp[j - padd_size] = 0;
pTmp[src_w + j] = 0;
} else { // reflect padding.
pTmp[j - padd_size] = pTmp[padd_size - j];
pTmp[src_w + padd_size - j] = pTmp[src_w - padd_size + j];
}
}
}
/* rows padding. */
if (padd_type == 0) {
memset(pBlur, 0, padd_size * img_w);
memset(pBlur + img_w * (img_h - 1 - padd_size), 0, padd_size * img_w);
} else {
for (int i = 0; i < padd_size; i++) {
memcpy(pBlur + i * img_w, pBlur + (k_size - 1 - i) * img_w, img_w);
memcpy(pBlur + (img_h - 1 - i) * img_w, pBlur + (img_h - k_size + i) * img_w, img_w);
}
}
/* Generate filter kernel. */
/* Eq: 1 / (2 * pi * Sigma^2) * exp( (x - m)^2 + (y - m)^2) */
kernelGenerate(gaussKernel, k_size, k_sigma);
/* Gaussian filter. */
for (int row = padd_size; row < img_h - padd_size; row++) {
for (int col = padd_size; col < img_w - padd_size; col++) {
float sum = 0;
for (int i = -padd_size; i < padd_size; i++) {
for (int j = -padd_size; j < padd_size; j++) {
sum += pBlur[(row + i) * img_w + col + j] * gaussKernel[(i + padd_size) * k_size + j + padd_size];
}
}
pDst[(row - padd_size) * src_w + col - padd_size] = sum;
}
}
free(gaussKernel);
free(pBlur);
}