【频率域平滑、锐化滤波器】理想滤波器，巴特沃思滤波器，高斯滤波器

一、背景知识
二、理想滤波器原理及实现
- 1.理想低通滤波器
- 2.理想低通滤波器的实现：
- 3.理想高通滤波器：
三、巴特沃思滤波器原理及实现
- 1.巴特沃思低通滤波器
- 2.巴特沃思高通滤波器
三、高斯滤波器原理及实现
- 1.高斯低通滤波器：
- 2.高斯高通滤波器：
四、代码附录
五、结尾

一、背景知识

本文主要介绍频率域滤波器，此处的频率域是基于傅立叶变换得出。

在一幅图像中，低频对应图像变化缓慢的部分，即图像大致外观和轮廓。高频部分对应图像变化剧烈的部分即图像细节。低通滤波器的功能是让低频率通过而滤掉或衰减高频，其作用是过滤掉包含在高频中的噪声。即低通滤波的效果是图像去噪声平滑增强，但同时也抑制了图像的边界即过滤掉图像细节，造成图像不同程序上的模糊。

这里考虑三种滤波器，分别为理想滤波器、巴特沃思滤波器和高斯滤波器。这三种滤波器分别有其低通滤波器和高通滤波器。

这三种滤波器涵盖了从非常尖锐（理想）到非常平坦（高斯）范围的滤波器函数。巴特沃思滤波器有一个参数，称为滤波器的“阶数”。当此参数的值较高时，巴特沃思滤波器接近理想滤波器。因此，巴特沃思滤波器可看作两种“极端”滤波器的过渡。

二、理想滤波器原理及实现

1.理想低通滤波器

我们所想像的最简单的低通滤波器是“截断”傅里叶变换中所有高频部分，这些成分处在距变换原点的距离比指定距离D0要更远得多的位置。这种滤波器被称为二维理想低通滤波器（ILPF），其变换函数为：

其中，D0是指定的非负数值，D（u,w）是（u，v）点距频率矩形原点的距离。如果要研究的图像尺寸为M*X，它的变换也有相同的尺寸，由于变换被中心化了，所以频率矩阵的中心在（u，v）=（M/2,N/2）处。在这种情况下，从点（u，v）到傅立叶变换中心的距离如下所示：（参照傅立叶变换和频率域）

需要注意理想低通滤波非常不实用，在使用中会产生“振铃”（输出图像的灰度剧烈变化处产生的震荡，就好像钟被敲击后产生的空气震荡）。然而，因为理想滤波器可以在计算机上实现，所以作为滤波概念发展的一部分，研究滤波器的特性非常有用。

2.理想低通滤波器的实现：

步骤一：将图像扩充为适合离散傅立叶变换的图像尺寸大小。

Mat Filter::ImageChangeForDft(const Mat & input)
{int width = getOptimalDFTSize(input.cols);int height = getOptimalDFTSize(input.rows);Mat output;copyMakeBorder(input,output,0,height-input.rows,0,width-input.cols,BORDER_CONSTANT,Scalar::all(0));output.convertTo(output,CV_32FC1);return output;
}

核心函数介绍：

int getOptimalDFTSize(int vecsize);

int类型的vecsize，向量尺寸，一般是图片的宽或高。

该函数可以自动计算最适合进行离散傅里叶变换的图像尺寸大小，从而提高计算速度。

CV_EXPORTS_W void copyMakeBorder(InputArray src, OutputArray dst, int top, int bottom, int left, int right,int borderType, const Scalar& value = Scalar() );

src：输入图像。 dst：输出图像(类型和输入图像类型一致，大小为：Size(src.cols+left+right,src.rows+top+bottom))。
top、bottom、left、 right: 这4个参数指定输出图像4个方向要扩展多少个像素，
bordertype：已拷贝的原图像长方形的边界的类型：

IPL_BORDER_CONSTANT：填充边界为固定值，值由函数最后一个参数指定。
IPL_BORDER_REPLICATE：边界用上下行或者左右列来复制填充。
value：如果border==BORDER_CONSTANT，value为border value。

步骤二：根据公式截断，获得滤波器ILPF

Mat Filter::IdealLowKernel(Mat &input,const float sigma) {Mat result(input.size(),CV_32FC1);float d0=sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小for(int i=0;i<input.rows ; i++ ){for(int j=0; j<input.cols ; j++ ){float d = sqrt(pow((i - input.rows/2),2) + pow((j - input.cols/2),2));//分子,计算pow必须为float型//d实际是到中心点的距离。if (d <= d0){result.at<float>(i,j)=1;}else{result.at<float>(i,j)=0;}}}return result; }

float d = sqrt(pow((i - input.rows/2),2) + pow((j - input.cols/2),2))

注意：该公式是基于频谱居中后的频谱图进行到中心点的距离计算，所以dft之后需要将频谱居中再计算

理想低通滤波“截断”傅里叶变换中所有高频部分，这些成分处在距变换原点的距离比指定距离D0要更远得多的位置，根据点（u，v）到傅立叶变换中心的距离进行截断。

步骤三：空间域转频率域，与滤波器ILPF相乘，再转回空间域

Mat Filter::FrequencyFilter(Mat &input,Mat &blur)
{//first, from spatial domain to frequency domain Mat plane[] = {input , Mat::zeros(input.size() , CV_32FC1)};Mat dftMat;merge(plane,2,dftMat);dft(dftMat,dftMat);//second, shift frequency domainsplit(dftMat,plane);
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0//频谱居中，四个象限位置互换，其中Plane[0]为实部，Plane[1]为虚部int cx = plane[0].cols/2;int cy = plane[0].rows/2;Mat part1_r(plane[0],Rect(0,0,cx,cy));Mat part2_r(plane[0],Rect(cx,0,cx,cy));Mat part3_r(plane[0],Rect(0,cy,cx,cy));Mat part4_r(plane[0],Rect(cx,cy,cx,cy));Mat temp;part1_r.copyTo(temp); part4_r.copyTo(part1_r);temp.copyTo(part4_r);part2_r.copyTo(temp); part3_r.copyTo(part2_r);temp.copyTo(part3_r);Mat part1_i(plane[1],Rect(0,0,cx,cy));  Mat part2_i(plane[1],Rect(cx,0,cx,cy));Mat part3_i(plane[1],Rect(0,cy,cx,cy));Mat part4_i(plane[1],Rect(cx,cy,cx,cy));part1_i.copyTo(temp); part4_i.copyTo(part1_i);temp.copyTo(part4_i);part2_i.copyTo(temp);  part3_i.copyTo(part2_i);temp.copyTo(part3_i);//third, multipy filter kernel and dftMatMat blur_r,blur_i,BLUR;multiply(plane[0], blur, blur_r); //滤波（实部与滤波器模板对应元素相乘）multiply(plane[1], blur,blur_i);//滤波（虚部与滤波器模板对应元素相乘）Mat plane1[]={blur_r, blur_i};merge(plane1,2,BLUR);//实部与虚部合并//last, from frequency domain to spatial domainidft( BLUR, BLUR);   split(BLUR,plane);magnitude(plane[0],plane[1],plane[0]); normalize(plane[0],plane[0],1,0,CV_MINMAX);Mat result;normalize(plane[0], result, 0, 255, CV_MINMAX);return result;
}

基于以上算法，可以更换滤波器达到不同的图像增强效果。

3.理想高通滤波器：

剔除低频，即当d<=d0 时为0。

Mat Filter::IdealHighKernel(Mat &scr,float sigma)
{Mat result(scr.size(),CV_32FC1);float d0=sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));//分子,计算pow必须为float型if (d <= d0){result.at<float>(i,j)=0;}else{result.at<float>(i,j)=1;}}}return result;
}

三、巴特沃思滤波器原理及实现

1.巴特沃思低通滤波器

n级巴特沃思低通滤波器（BLPF）的传递函数（且截止频率距原点的距离为D0）的定义如下：

不同于ILPF，BLPF变换函数在通带与被滤除的频率之间没有明显的截断。对于有传递函数的滤波器，定义一个截止频率的位置并使H(u,v)幅度降到其最大值的一部分。在式中，当D(u,v)=D0时，H(u,v)=0.5（从最大值1降到它的50%）

一个一阶的巴特沃思滤波器没有振铃。在二阶中振铃通常很微小，但阶数增高的时候振铃变会成为一个重要的影响因素。阶数为5明显振铃和负值，阶数为20与ILPF相似。二阶的BLPF是在有效的低通滤波和可接受的振铃特性之间的折中。

Mat Filter::ButterworthLowKernel(Mat &scr,float sigma, int n)
{Mat result(scr.size(),CV_32FC1);double D0 = sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));result.at<float>(i,j)=1.0 / (1 + pow(d / D0, 2 * n));}}string name = "ButterworthLowKerneld0=" + std::to_string(sigma) + "n=" + std::to_string(n);imshow(name, result);return result;
}

2.巴特沃思高通滤波器

同低通滤波器情况一样，可以认为巴特沃思高通滤波器比IHPF更平滑。

Mat Filter::ButterworthHighKernel(Mat &scr,float sigma, int n)
{Mat result(scr.size(),CV_32FC1); //，CV_32FC1double D0 = sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));result.at<float>(i,j)=1.0 / (1 + pow(D0 / d, 2 * n));}}return result;
}

三、高斯滤波器原理及实现

1.高斯低通滤波器：

其中，D0是截止频率，D(u,v)是距频率矩形中心的距离。高斯滤波器的宽度由参数 D0 表征，决定了平滑程度，而且 D0越大，高斯滤波器的频带就越宽，平滑程度就越好。因为噪声主要集中在高频段，所以通过高斯低通滤波器可以滤除噪声信息、平滑图像，但与此同时会滤除图像的细节信息，使图像变得模糊。

Mat Filter::GaussianLowPassKernel(Mat scr,float sigma)
{Mat gaussianBlur(scr.size(),CV_32FC1);float d0=2*sigma*sigma;//高斯函数参数，越小，频率高斯滤波器越窄，滤除高频成分越多，图像就越平滑for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){float d=pow(float(i-scr.rows/2),2)+pow(float(j-scr.cols/2),2);//分子,计算pow必须为float型gaussianBlur.at<float>(i,j)=expf(-d/d0);//expf为以e为底求幂（必须为float型）}}return gaussianBlur;
}

2.高斯高通滤波器：

通过高斯高通滤波器可以增强细节信息，提升图像的高频分量，减少低频分量，对微小物体和细线条也能很好地增强显示。

Mat Filter::GaussianHighPassKernel(Mat scr,float sigma)
{Mat gaussianBlur(scr.size(),CV_32FC1); //，CV_32FC1float d0=2*sigma*sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){float d=pow(float(i-scr.rows/2),2)+pow(float(j-scr.cols/2),2);//分子,计算pow必须为float型gaussianBlur.at<float>(i,j)=1-expf(-d/d0);}}return gaussianBlur;
}

通过公式里的傅立叶反变换得到的空间高斯滤波器将没有振铃。也是平时最为常用的滤波器。

四、代码附录

<filter.h>

#ifndef FILTER
#define FILTER#include <iostream>
#include <opencv2/opencv.hpp>class Filter
{public:cv::Mat ImageChangeForDft(const cv::Mat & input);cv::Mat IdealLowKernel(cv::Mat &input,const float sigma);cv::Mat FrequencyFilter(cv::Mat &input,cv::Mat &blur);cv::Mat IdealLowPassFilter(const cv::Mat &input, float sigma);cv::Mat ButterworthLowKernel(cv::Mat &scr,float sigma, int n);cv::Mat ButterworthLowPassFilter(cv::Mat &src, float d0, int n);cv::Mat GaussianLowPassKernel(cv::Mat scr,float sigma);cv::Mat GaussianLowPassFilter(cv::Mat &src, float d0);cv::Mat IdealHighKernel(cv::Mat &scr,float sigma);cv::Mat IdealHighPassFilter(cv::Mat &src, float sigma);cv::Mat ButterworthHighKernel(cv::Mat &scr,float sigma, int n);cv::Mat ButterworthHighPaassFilter(cv::Mat &src, float d0, int n);cv::Mat GaussianHighPassKernel(cv::Mat scr,float sigma);cv::Mat GaussianHighPassFilter(cv::Mat &src, float d0);};

<filter.cpp>

#include <iostream>
#include <opencv2/opencv.hpp>#include "filter.h"using namespace std;
using namespace cv;Mat Filter::ImageChangeForDft(const Mat & input)
{int width = getOptimalDFTSize(input.cols);int height = getOptimalDFTSize(input.rows);Mat output;copyMakeBorder(input,output,0,height-input.rows,0,width-input.cols,BORDER_CONSTANT,Scalar::all(0));output.convertTo(output,CV_32FC1);return output;
}//*****************理想低通滤波器******************
Mat Filter::IdealLowKernel(Mat &input,const float sigma)
{Mat result(input.size(),CV_32FC1);float d0=sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小for(int i=0;i<input.rows ; i++ ){for(int j=0; j<input.cols ; j++ ){float d = sqrt(pow((i - input.rows/2),2) + pow((j - input.cols/2),2));//分子,计算pow必须为float型//d实际是到中心点的距离。if (d <= d0){result.at<float>(i,j)=1;}else{result.at<float>(i,j)=0;}}}return result;
}Mat Filter::FrequencyFilter(Mat &input,Mat &blur)
{//first, from spatial domain to frequency domain Mat plane[] = {input , Mat::zeros(input.size() , CV_32FC1)};Mat dftMat;merge(plane,2,dftMat);dft(dftMat,dftMat);//second, shift frequency domainsplit(dftMat,plane);
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0//频谱居中，四个象限位置互换，其中Plane[0]为实部，Plane[1]为虚部int cx = plane[0].cols/2;int cy = plane[0].rows/2;Mat part1_r(plane[0],Rect(0,0,cx,cy));Mat part2_r(plane[0],Rect(cx,0,cx,cy));Mat part3_r(plane[0],Rect(0,cy,cx,cy));Mat part4_r(plane[0],Rect(cx,cy,cx,cy));Mat temp;part1_r.copyTo(temp); part4_r.copyTo(part1_r);temp.copyTo(part4_r);part2_r.copyTo(temp); part3_r.copyTo(part2_r);temp.copyTo(part3_r);Mat part1_i(plane[1],Rect(0,0,cx,cy));  Mat part2_i(plane[1],Rect(cx,0,cx,cy));Mat part3_i(plane[1],Rect(0,cy,cx,cy));Mat part4_i(plane[1],Rect(cx,cy,cx,cy));part1_i.copyTo(temp); part4_i.copyTo(part1_i);temp.copyTo(part4_i);part2_i.copyTo(temp);  part3_i.copyTo(part2_i);temp.copyTo(part3_i);//third, multipy filter kernel and dftMatMat blur_r,blur_i,BLUR;multiply(plane[0], blur, blur_r); //滤波（实部与滤波器模板对应元素相乘）multiply(plane[1], blur,blur_i);//滤波（虚部与滤波器模板对应元素相乘）Mat plane1[]={blur_r, blur_i};merge(plane1,2,BLUR);//实部与虚部合并//last, from frequency domain to spatial domainidft( BLUR, BLUR);   split(BLUR,plane);magnitude(plane[0],plane[1],plane[0]); normalize(plane[0],plane[0],1,0,CV_MINMAX);Mat result;normalize(plane[0], result, 0, 255, CV_MINMAX);return result;
}//理想低通滤波器
Mat Filter::IdealLowPassFilter(const Mat &input, float sigma)
{Mat padded = ImageChangeForDft(input);Mat ideal_kernel=IdealLowKernel(padded,sigma);Mat result = FrequencyFilter(padded,ideal_kernel);return result;
}Mat Filter::ButterworthLowKernel(Mat &scr,float sigma, int n)
{Mat result(scr.size(),CV_32FC1);double D0 = sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));result.at<float>(i,j)=1.0 / (1 + pow(d / D0, 2 * n));}}string name = "ButterworthLowKerneld0=" + std::to_string(sigma) + "n=" + std::to_string(n);imshow(name, result);return result;
}//巴特沃斯低通滤波器
Mat Filter::ButterworthLowPassFilter(Mat &src, float d0, int n)
{//H = 1 / (1+(D/D0)^2n)    n表示巴特沃斯滤波器的次数//阶数n=1 无振铃和负值    阶数n=2 轻微振铃和负值  阶数n=5 明显振铃和负值   阶数n=20 与ILPF相似Mat padded = ImageChangeForDft(src);Mat butterworth_kernel=ButterworthLowKernel(padded,d0, n);Mat result = FrequencyFilter(padded,butterworth_kernel);return result;
}Mat Filter::GaussianLowPassKernel(Mat scr,float sigma)
{Mat gaussianBlur(scr.size(),CV_32FC1);float d0=2*sigma*sigma;//高斯函数参数，越小，频率高斯滤波器越窄，滤除高频成分越多，图像就越平滑for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){float d=pow(float(i-scr.rows/2),2)+pow(float(j-scr.cols/2),2);//分子,计算pow必须为float型gaussianBlur.at<float>(i,j)=expf(-d/d0);//expf为以e为底求幂（必须为float型）}}return gaussianBlur;
}//高斯低通
Mat Filter::GaussianLowPassFilter(Mat &src, float d0)
{Mat padded = ImageChangeForDft(src);Mat gaussian_kernel=GaussianLowPassKernel(padded,d0);//理想低通滤波器Mat result = FrequencyFilter(padded,gaussian_kernel);return result;
}Mat Filter::IdealHighKernel(Mat &scr,float sigma)
{Mat result(scr.size(),CV_32FC1);float d0=sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));//分子,计算pow必须为float型if (d <= d0){result.at<float>(i,j)=0;}else{result.at<float>(i,j)=1;}}}return result;
}//理想高通滤波器
Mat Filter::IdealHighPassFilter(Mat &src, float sigma)
{Mat padded = ImageChangeForDft(src);Mat ideal_kernel=IdealHighKernel(padded,sigma);Mat result = FrequencyFilter(padded,ideal_kernel);return result;
}Mat Filter::ButterworthHighKernel(Mat &scr,float sigma, int n)
{Mat result(scr.size(),CV_32FC1); //，CV_32FC1double D0 = sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){double d = sqrt(pow((i - scr.rows/2),2) + pow((j - scr.cols/2),2));result.at<float>(i,j)=1.0 / (1 + pow(D0 / d, 2 * n));}}return result;
}//巴特沃斯高通滤波器
Mat Filter::ButterworthHighPaassFilter(Mat &src, float d0, int n)
{//H = 1 / (1+(D0/D)^2n)    n表示巴特沃斯滤波器的次数Mat padded = ImageChangeForDft(src);Mat butterworth_kernel=ButterworthHighKernel(padded,d0, n);Mat result = FrequencyFilter(padded,butterworth_kernel);return result;
}Mat Filter::GaussianHighPassKernel(Mat scr,float sigma)
{Mat gaussianBlur(scr.size(),CV_32FC1); //，CV_32FC1float d0=2*sigma*sigma;for(int i=0;i<scr.rows ; i++ ){for(int j=0; j<scr.cols ; j++ ){float d=pow(float(i-scr.rows/2),2)+pow(float(j-scr.cols/2),2);//分子,计算pow必须为float型gaussianBlur.at<float>(i,j)=1-expf(-d/d0);}}return gaussianBlur;
}//高斯高通
Mat Filter::GaussianHighPassFilter(Mat &src, float d0)
{Mat padded = ImageChangeForDft(src);Mat gaussian_kernel=GaussianHighPassKernel(padded,d0);//理想低通滤波器Mat result = FrequencyFilter(padded,gaussian_kernel);return result;
}#endif

五、结尾

至此还可以引出其他拓展点：
如频率域的锐化滤波器还有拉普拉斯算子等、同时还有钝化模板、高频提升滤波和高频加强滤波、同态滤波器等

【图像处理：频率域平滑与锐化】理想滤波器，巴特沃思滤波器，高斯滤波器相关推荐

数字图像处理频率域平滑 MATLAB实验
一.原理_频率域平滑理想低通滤波器的传递函数为: n 阶巴特沃斯低通滤波器的传递函数为: n 阶指数低通滤波器的传递函数为: 二.步骤 (1)读入原图像test.tif并显示: (2)对原图像添加高 ...
数字图像处理-频率域滤波
一.实验原理频率域滤波是对图像进行傅里叶变换,将图像由图像空间转换到频域空间,然后在频率域中对图像的频谱作分析处理,以改变图像的频率特征,原理是用傅里叶变换表示的函数特征完全可以通过傅里叶反变换来重 ...
数字图像处理频率域锐化 MATLAB实验
一.原理_频率域锐化理想高通滤波器的传递函数为: n 阶巴特沃斯高通滤波器的传递函数为: n 阶指数高通滤波器的传递函数为: 二.步骤 (1)读入原图像test.tif并显示: (2)采用理想高通滤 ...
计算机图像处理实验三图像空间域平滑与锐化
一.实验目的与要求 1．加深对图像增强及边缘检测技术的感性认识,应用MATLAB工具箱自带的处理函数或自己编程完成相关的工作,分析处理结果,巩固所学理论知识. 2．熟练掌握空域滤波中常用的平滑和锐化滤 ...
数字图像处理-频率域滤波原理
from:https://blog.csdn.net/forrest02/article/details/55510711?locationNum=15&fps=1 写在前面的话作者是一名在 ...
图像处理-空间域平滑滤波
个人博客:http://www.chenjianqu.com/ 原文链接:http://www.chenjianqu.com/show-12.html 目录: 1. 空间域和频域的概念 2 ...
数字图像处理空间域平滑 MATLAB实验
一.原理_均值滤波设加噪图像为 f(x,y) ,经均值滤波处理后的图像为g(x,y) ,则: 式中, S是(x,y)像素点的邻域,K是S内的像素数. 最典型的S为3X3邻域,可用模板形式表示为: 这 ...
python实现陷波滤波器、低通滤波器、高斯滤波器、巴特沃斯滤波器
在一幅图像中,其低频成分对应者图像变化缓慢的部分,对应着图像大致的相貌和轮廓,而其高频成分则对应着图像变化剧烈的部分,对应着图像的细节(图像的噪声也属于高频成分). 滤波器低通滤波器高通滤波器陷 ...
数字图像处理——第四章频率域图像增强
文章目录频率域图像增强 1. 傅里叶变换 1.1 一维傅里叶变换 1.2 二维傅里叶变换 2. 频率域滤波 2.1 陷波滤波器及其性质 2.2 空间域滤波与频率域滤波之间的对应关系 2.2.1 空间 ...

【图像处理：频率域平滑与锐化】理想滤波器，巴特沃思滤波器，高斯滤波器