0. 完整仓库 & 教程:
这一章节的完整代码在:Chapter 7. Image Noise and Noise Reduction
如果你喜欢这个系列的文章或者感觉对你有帮助,请给我的仓库一个⭐️。
1. 同态滤波器
算法:
同态滤波器的过程与第六章(这部分可以参考:6:低通滤波器和高通滤波器)中的高通/低通滤波器类似:
首先对图像进行二维离散傅里叶变换 DFT (这部分原理需要理解: 5: 二维离散傅里叶变换和重建),得到每个像素的复数数组。 每个复数的实部和虚部分别乘以 𝐻(𝑢, 𝑣)。 𝐻(𝑢, 𝑣) 的公式为:
其中 𝛾𝐻 > 1, 𝛾𝐿 < 1,并且 表示通带的半径。 𝐷(𝑢, 𝑣)的计算方法也是两点之间的距离,通过公式:
其中 𝑃 和𝑄 是图像的长和宽。在后续的代码处理中,设置𝛾𝐻 = 1.5, 𝛾𝐿 = 0.75, = 30, 𝑐 = 1。
图像和结果对比(bridge, goldhill):
结果分析:
同态滤波本质上是频域滤波器(frequency domain filter),它结合了频率滤波和空间灰度变换。 通常,图像的亮部(illumination)是缓慢变化的部分,集中在低频部分。 反射部分(reflection)属于快速变化的部分,往往集中在高频部分。
这里选用模型的形状类似于一个逆高斯滤波器(见章节6),可以减弱低频分量(illumination)的影响并强调高频分量(reflection)。 最终的结果是压缩了动态范围和增强了对比度,使纹理更加突出。
代码实现(完整代码见顶部GitHub):
for(int i = 0; i < height; i++) {
for(int j = 0; j < width; j++) {
float dis = sqrt(pow((float)i - height/2, 2) + pow((float)j - width/2, 2));
float H = (1.5 - 0.75) * (1 - pow(M_E, pow(dis / 30, 2) * -1)) + 0.75;
real[(int)(i*height + j)] = real_array[(int)(i*height + j)] * H;
imaginary[(int)(i*height + j)] = imaginary_array[(int)(i*height + j)] * H;
}
}
printf("Homomorphic Finished!\n");
outimage = iDFT(image, real, imaginary);
这里再把第5章里傅里叶变换和重建的代码搬过来下:
void DFT(Image *image, float *real_array, float *imaginary_array) {
unsigned char *tempin, *tempout;
float real, imaginary;
Image *outimage;
outimage = CreateNewImage(image, (char*)"#testing function");
tempin = image->data;
tempout = outimage->data;
printf("DFT algorithm is executing...\n");
for(int i = 0; i < image->Height; i++) {
for(int j = 0; j < image->Width; j++) {
real = 0;
imaginary = 0;
for(int m = 0; m < image->Height; m++) {
for(int n = 0; n < image->Width; n++) {
float temp = (float)i * m / (float)image->Height + (float)j * n / (float)image->Width;
int offset = (m + n) % 2 == 0 ? 1 : -1;
real += tempin[image->Width * m + n] * cos(-2 * pi * temp) * offset;
imaginary += tempin[image->Width * m + n] * sin(-2 * pi * temp) * offset;
}
}
real_array[image->Width * i + j] = real;
imaginary_array[image->Width * i + j] = imaginary;
int temp = (int)(sqrt(real*real + imaginary*imaginary) / sqrt((float)image->Height*(float)image->Width));
if(temp < 0) temp = 0;
if(temp > 255) temp = 255;
tempout[image->Width * i + j] = temp;
}
}
printf("DFT Finished!\n");
SavePNMImage(outimage, (char*)"DFT.pgm");
}
Image *iDFT(Image *image, float *real_array, float *imaginary_array) {
unsigned char *tempin, *tempout;
float real;
Image *outimage;
outimage = CreateNewImage(image, (char*)"#testing fun