文章目录

1. 频域中图像平滑流程
2. 理想低通滤波器
3. 巴特沃斯低通滤波器
4. 指数低通滤波器
5. 梯形低通滤波器

1. 频域中图像平滑流程

在一副图像中，平滑的图像信号在频谱中贡献低频分量，而高频部分往往是图像的细节和边界，噪声也具有丰富的高频分量。所以对图像进行低通滤波可以消减噪声，但也可能使边界和图像细节变得模糊。
频域中图像滤波的处理流层如下图：

其中DFT为离散傅里叶变换，IDFT为离散傅里叶反变换，H(m,n)H(m,n)H(m,n)为滤波器的传递函数。滤波器的输出为：

2. 理想低通滤波器

二维理想低通滤波器(Ideal Low Pass Filter,ILPF )的传递函数定义为：（连续形式给出）

其中D(u,v)=(u2+v2)12，D0D(u,v)=(u^2+v^2)^\frac{1}{2}，D0D(u,v)=(u2+v2)21，D0为截止频率。H(u,v)H(u,v)H(u,v)的三维示意图如下：

二维剖面如下：

传递函数H(u,v)H(u,v)H(u,v)表示以D0D0D0为半径的圆域内的分量无损的通过，其他分量则完全的被滤除掉。
对于数字图像，D(u,v)D(u,v)D(u,v)可以用离散的D(m,n)D(m,n)D(m,n)代替，D(m,n)D(m,n)D(m,n)同样表示为某频率分量离零频率分量的距离。因为图像进行DFT后零频率分量在频谱图的四个角，为了观察和操作方便我们需要把零频率分量移到频谱中心，操作完后再移回原位。
下图是一副数字图像的低通滤波示意图：（零频率分量已经移到频谱中心，中心最亮的点既是零频率分量）

对其进行低通滤波，相当于圆内的频率通过，圆外滤除。
理想低通滤波器关键代码：

Mat ImageSmoothing::LowFilter(Mat input, double D0)//理想低通滤波
{//D0为低通滤波器的截止频率if (input.channels() != 1) {std::cout << "LowFilter只能处理单通道图像，要处理多通道图像请用:LowFilter3D" << std::endl;exit(0);}Mat src = input.clone();Mat srcDft = DFT::dft2(src);//傅里叶变换srcDft = DFT::shift(srcDft);//将零频率分量移到频谱中心int M = srcDft.rows;int N = srcDft.cols;int m1 = M / 2;int n1 = N / 2;int channels = srcDft.channels();for (int m = 0; m < M; m++) {double* current = srcDft.ptr<double>(m);for (int n = 0; n < N; n++) {if (sqrt((m - m1) * (m - m1) + (n - n1) * (n - n1)) > D0) {current[n*channels] = 0.0;current[n * channels+1] = 0.0;}}}srcDft = DFT::shift(srcDft);//还原Mat dst;idft(srcDft, dst, DFT_SCALE | DFT_REAL_OUTPUT);dst.convertTo(dst, CV_8U);return dst;
}

其中 DFT::shift与 DFT::dft2定义如下：

Mat DFT::shift(Mat input)
{Mat src = input.clone();src = src(Range(0, src.rows & -2), Range(0, src.cols & -2));int y0 = src.rows / 2;int x0 = src.cols / 2;Mat q0(src, Rect(0, 0, x0, y0));       //左上角Mat q1(src, Rect(x0, 0, x0, y0));      //右上角Mat q2(src, Rect(0, y0, x0, y0));      //左下角Mat q3(src, Rect(x0, y0, x0, y0));     //右下角swapArea(q0, q3);swapArea(q1, q2);return src;
}
Mat DFT::dft2(Mat input)
{if (input.channels() != 1) {std::cout << "只处理单通道的图像" << std::endl;exit(0);}Mat src = input.clone();//扩充图像尺寸为DFT的最佳尺寸Mat padded;int newRows = getOptimalDFTSize(src.rows);int newCols = getOptimalDFTSize(src.cols);copyMakeBorder(src, padded, newRows - src.rows, 0, newCols - src.cols, 0, BORDER_CONSTANT, Scalar(0));//将planes融合合并成一个多通道数组Mat plans[] = { Mat_<double>(padded),Mat::zeros(padded.size(),CV_64F) };Mat mergeArray;merge(plans, 2, mergeArray);//傅里叶变换dft(mergeArray, mergeArray);return mergeArray;
}

上面的代码只能实现单通道图像，对于多通道图像，可以使用上面的LowFilter函数对其每一个通道分别处理即可，代码如下：

Mat ImageSmoothing::LowFilter3D(Mat input, double D0)
{Mat src = input.clone();int channels = src.channels();Mat *dst=new Mat[channels];split(src, dst);for (int i = 0; i < channels; i++) {dst[i] = LowFilter(dst[i],D0);}Mat dstImg;merge(dst, channels, dstImg);return dstImg;
}

理想低通滤波效果：（取D0=100D0=100D0=100）

理想低通滤波器由于是锐截止，处理后的图像会出现不应有的亮环，也就是“振铃效应”，而且图像也会变得模糊一些。

3. 巴特沃斯低通滤波器

巴特沃斯低通滤波器（Butterworth Low Pass Filter,BLPF）的传递函数如下：

其中D(u,v)=(u2+v2)12D(u,v)=(u^2+v^2)^\frac{1}{2}D(u,v)=(u2+v2)21，D0D_0D0为截止频率，nnn为阶数，取正整数，控制衰减速度。H(u,v)H(u,v)H(u,v)的三维示意图与二维剖面图如下：

这种滤波器不是锐截止，无“振铃”现象，提高了图像的细节清晰度。
巴特沃斯低通滤波器关键代码：

Mat ImageSmoothing::ButterworthLPF(Mat input, double D0, int n0)
{//D0为巴特沃斯滤波器的截止频率//n0控制衰减速度if (input.channels() != 1) {std::cout << "ButterworthLPF只能处理单通道图像，要处理多通道图像请用:ButterworthLPF3D" << std::endl;exit(0);}Mat src = input.clone();Mat srcDft = DFT::dft2(src);srcDft = DFT::shift(srcDft);//将零频率分量移到频谱中心int M = srcDft.rows;int N = srcDft.cols;int m1 = M / 2;int n1 = N / 2;int channels = srcDft.channels();for (int m = 0; m < M; m++) {double* current = srcDft.ptr<double>(m);for (int n = 0; n < N; n++) {double D= sqrt((m - m1) * (m - m1) + (n - n1) * (n - n1));double H = 1 / (1 + (sqrt(2) - 1) *pow(D / D0,2*n0));current[n * channels] *= H;current[n * channels + 1] *= H;}}srcDft = DFT::shift(srcDft);//还原Mat dst;idft(srcDft, dst, DFT_SCALE | DFT_REAL_OUTPUT);dst.convertTo(dst, CV_8U);return dst;
}

同样的，对于多通道图像，可以使用上面的ButterworthLPF函数对其每一个通道分别处理即可，代码如下：

Mat ImageSmoothing::ButterworthLPF3D(Mat input, double D0, int n0)
{Mat src = input.clone();int channels = src.channels();Mat* dst = new Mat[channels];split(src, dst);for (int i = 0; i < channels; i++) {dst[i] = ButterworthLPF(dst[i], D0,n0);}Mat dstImg;merge(dst, channels, dstImg);return dstImg;
}

巴特沃斯低通滤波效果：（取D0=100，n=2D0=100，n=2D0=100，n=2）

4. 指数低通滤波器

指数滤波器的传递函数为：

其中D(u,v)=(u2+v2)12D(u,v)=(u^2+v^2)^\frac{1}{2}D(u,v)=(u2+v2)21，D0D_0D0为截止频率，nnn为阶数，取正整数，控制衰减速度。H(u,v)H(u,v)H(u,v)的三维示意图与二维剖面图如下：

指数低通滤波器关键代码：

Mat ImageSmoothing::expLPF(Mat input, double D0, int n0)
{//指数低通滤波器if (input.channels() != 1) {std::cout << "expLPF只能处理单通道图像，要处理多通道图像请用:expLPF3D" << std::endl;exit(0);}Mat src = input.clone();Mat srcDft = DFT::dft2(src);srcDft = DFT::shift(srcDft);//将零频率分量移到频谱中心int M = srcDft.rows;int N = srcDft.cols;int m1 = M / 2;int n1 = N / 2;int channels = srcDft.channels();for (int m = 0; m < M; m++) {double* current = srcDft.ptr<double>(m);for (int n = 0; n < N; n++) {double D = sqrt((m - m1) * (m - m1) + (n - n1) * (n - n1));double H = exp((log(1/sqrt(2)))*pow(-D/D0,n0));current[n * channels] *= H;current[n * channels + 1] *= H;}}srcDft = DFT::shift(srcDft);//还原Mat dst;idft(srcDft, dst, DFT_SCALE | DFT_REAL_OUTPUT);dst.convertTo(dst, CV_8U);return dst;
}

同样的，对于多通道图像，可以使用上面的expLPF函数对其每一个通道分别处理即可，代码如下：

Mat ImageSmoothing::expLPF3D(Mat input, double D0, int n0)
{Mat src = input.clone();int channels = src.channels();Mat* dst = new Mat[channels];split(src, dst);for (int i = 0; i < channels; i++) {dst[i] = expLPF(dst[i], D0, n0);}Mat dstImg;merge(dst, channels, dstImg);return dstImg;
}

指数低通滤波效果：（取D0=100，n=2D0=100，n=2D0=100，n=2）

5. 梯形低通滤波器

梯形滤波器的传递函数为：

其中D(u,v)=(u2+v2)12D(u,v)=(u^2+v^2)^\frac{1}{2}D(u,v)=(u2+v2)21，D0D_0D0为截止频率，D0D_0D0和D1D_1D1按图像给定的特点预先给定，且D0<D1D_0<D_1D0<D1。H(u,v)H(u,v)H(u,v)的三维示意图与二维剖面图如下：

梯形低通滤波器关键代码：

Mat ImageSmoothing::trapezoidLPF(Mat input, double D0, double D1)
{//梯形低通滤波if (input.channels() != 1) {std::cout << "trapezoidLPF只能处理单通道图像，要处理多通道图像请用:trapezoidLPF3D" << std::endl;exit(0);}if (D0>=D1) {std::cout << "输入参数有误，参数必须是D0<D1" << std::endl;exit(0);}Mat src = input.clone();Mat srcDft = DFT::dft2(src);srcDft = DFT::shift(srcDft);//将零频率分量移到频谱中心int M = srcDft.rows;int N = srcDft.cols;int m1 = M / 2;int n1 = N / 2;int channels = srcDft.channels();for (int m = 0; m < M; m++) {double* current = srcDft.ptr<double>(m);for (int n = 0; n < N; n++) {double D = sqrt((m - m1) * (m - m1) + (n - n1) * (n - n1));double H;if (D < D0) {H = 1;}else if (D > D1) {H = 0;}else {H = (D - D1) / (D0 - D1);}current[n * channels] *= H;current[n * channels + 1] *= H;}}srcDft = DFT::shift(srcDft);//还原Mat dst;idft(srcDft, dst, DFT_SCALE | DFT_REAL_OUTPUT);dst.convertTo(dst, CV_8U);return dst;
}

同样的，对于多通道图像，可以使用trapezoidLPF函数对其每一个通道分别处理即可，代码如下：

Mat ImageSmoothing::trapezoidLPF3D(Mat input, double D0, double D1)
{Mat src = input.clone();int channels = src.channels();Mat* dst = new Mat[channels];split(src, dst);for (int i = 0; i < channels; i++) {dst[i] = trapezoidLPF(dst[i], D0, D1);}Mat dstImg;merge(dst, channels, dstImg);return dstImg;
}

梯形低通滤波效果：（取D0=100，D1=150D_0=100，D_1=150D0=100，D1=150）

频域中图像平滑的方法相关推荐

UA OPTI501 电磁波求解麦克斯韦方程组的Fourier方法1 在频域中讨论麦克斯韦方程组
UA OPTI501 电磁波求解麦克斯韦方程组的Fourier方法1 在频域中讨论麦克斯韦方程组标量势与矢量势 Lorenz Gauge 麦克斯韦方程组的意义是阐述了电磁场的源(ρfree,Jfr ...
UA OPTI501 电磁波经典电动力学中的Fourier方法基础
UA OPTI501 电磁波经典电动力学中的Fourier方法基础 4-D Spatial-temporal Fourier变换常用的Fourier变换结论 1的Fourier变换是Dirac函数 ...
【Python】这篇文章能让你明白经验模态分解（EMD）——EMD在python中的实现方法
暂时打断一下滤波专题,插播一条EMD在python中实现方法的文章. 本篇是Mr.看海:这篇文章能让你明白经验模态分解(EMD)--EMD在MATLAB中的实现方法的姊妹篇,也就是要在python中实 ...
射频原理和组件介绍：学习如何在频域中分析射频（RF）问题
什么是频域?为什么它对射频(RF)设计,分析和测试如此有价值? 熟练掌握射频(RF)设计过程中最基本的步骤之一就是学会在频域中进行问题思考.对于我们大多数人来说,我们在电路和信号方面的早期经验绝大部分 ...
matlab频域滤波处理，手动滤波，图像去网格背景：频域中手动选择目标区域进行滤波，滤掉噪声，灵活地达到理想的效果。
目录一.原理简介二.现有方法三.代码实现四.测试结果图一.原理简介主要思路为:通过傅里叶变换,将图像显示在频率域中,通过观察频谱图,如果了解频谱图与原图之间的关联,可以较为容易地发现噪声在 ...
信号与线性系统中LTI系统分析方法的不断优化
在信号与线性系统分析中,分析对象是线性时不变系统(LTI),即输入信号的线性组合等于输出信号的线性组合,且系统的参数不随时间变化.对应描述该系统的数学模型可以是常系数微分方程(输入信号为连续信号),常 ...
vue打印props的值_vue中props传值方法
vue中props传值方法 1.开发环境 vue 2.电脑系统 windows10专业版 3.在使用 vue开发的过程中,我们经常会使用 props进行组件的传值,下面是我的分享,希望对你有所帮助! ...
Javascript中二进制数据处理方法
Javascript中二进制数据处理方法转载于:https://www.cnblogs.com/motadou/archive/2012/02/19/2358514.html
HashSet中的add()方法( 五 )（详尽版）
上接 HashSet中的add()方法( 四 )(详尽版) ,我们再重写一下equals()方法来看看是否可以不能存入相同的id: 在学生类中再重写equals()方法: public class S ...

频域中图像平滑的方法