数字图像处理第四章频率域滤波（低通滤波器、高通滤波器、拉普拉斯滤波、同态滤波器）

最新推荐文章于 2023-08-23 00:15:00 发布
最新推荐文章于 2023-08-23 00:15:00 发布
阅读量1.3w 收藏 85
点赞数 9
分类专栏： 传统图像处理 文章标签： opencv 计算机视觉 图像处理
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/wangjia2575525474/article/details/116993230
版权

本章节的主要内容具体包括：傅里叶变换的概念及处理的相关知识、频率域卷积概念、三种低通滤波器的原理及代码实现、三种高通滤波器的原理及代码实现、频率域拉普拉斯算法原理及实现、同态滤波器原理及代码实现。

4.1傅里叶变换原理

频率域图像处理步骤：
在这里插入图片描述
在具体进行频率域的各种处理滤波的前后,进行了傅立叶变换以及傅立叶反变换.这两个变换的过程就是将空间的信息分解为在频率上的表示,或者将频率上的表示转化为空间上的表示,两种变换是互为逆变换的.正是通过傅立叶正反变换的处理,才使得频率域上的处理可以用于图像的增强。

1、傅里叶变换

概念：
法国数学家傅里叶指出：任何周期函数都可以表示为不同频率的正弦和或余弦和的形式，每个正弦和或者余弦和乘以不同的系数（这个正弦和或余弦和就是傅里叶级数）。
非周期函数也可以用正弦或者余弦乘以加权函数的积分来表示。这种情况下就是傅里叶变换。。DFT指的是单变量的离散傅里叶变化；
公式：（二维离散形式）
正变换：
在这里插入图片描述
反变换：
在这里插入图片描述
图像的频率是表征图像中灰度变化剧烈程度的指标，是灰度在平面空间上的梯度。直流分量 [公式] 对应一幅图像的平均灰度，低频对应图像中变化缓慢的灰度分量，高频部分对应图像边缘和灰度级突变的部分。

2、频率域卷积–卷积定理

首先在空间域中，两个函数的卷积涉及一个函数关于原点做翻转（旋转180度）并划过另一个函数，在滑动在过程中在每一个位移处我们都做计算，在空间域中是做乘积之和；
而在频率域中，两个连续函数在卷积变为求积分：
在这里插入图片描述
然后我们利用傅里叶变换将函数扩展到无穷大：
在这里插入图片描述
其中将t所在在域成为空间域，将u所在的予称为频率域
结论：空间域两个函数的卷积的傅里叶变换 = 两个函数傅里叶变换在频率域中在乘积
如果有两个函数傅里叶变换在乘积，那么我们可以通过计算傅里叶反变换而得到空间域卷积
f与F互为傅里叶变化对
相反的，如果有空间域中在两个函数在乘积我们可以通过傅里叶变换得到频率域卷积

3、图像频域中心化处理

用傅里叶变换的平移特性证明
在这里插入图片描述

其中
在这里插入图片描述在这里插入图片描述
由傅里叶变换的平移特性可得
在这里插入图片描述

因此 在这里插入图片描述

去中心化处d理的原因：一张图像中的能量主要集中在低频区，中心化处理可以将图像的低频分量集中在频谱图中央，便于观察。
中心化：用（-1^x）乘以f(x)将位于原点在数据F(0)移动到【0，M-1】的中心位置。在这里插入图片描述

4.2二维离散傅里叶变换的一些性质

1、平移和旋转

用指数项乘以f(x,y)将使得DFT的原点移到（u0,v0）;反之，用负指数乘以F（u,v）将使f(x,y)的原点移到点（x0,y0）
平移不影响F(U,V)的幅度（谱）
在这里插入图片描述

旋转：F（u,v）旋转一个角度，则f(x,y)旋转相同在角度，反之亦然。

2.周期性

二维傅里叶变换和反变换在U和v方向是无限周期的。
在这里插入图片描述
中心化：用（-1^x）乘以f(x)将位于原点在数据F(0)移动到【0，M-1】的中心位置。

在这里插入图片描述

3、对称性

4、频率域滤波

低频对应图像中变化缓慢的灰度分量，对应着图像总体灰度级的显示
高频部分对应图像边缘和灰度级突变的部分，对应着图像中变化快在分量，图像的细节。
频率域滤波的步骤：
在这里插入图片描述
对于不同的滤波器，滤波函数H(u,v)不相同；
滤波器有陷波滤波器：除了原点有凹陷外，其他均是常量函数；，此滤波器可以设置F(0，0)为0，而保留其他傅里叶变换在频率成分不变。处理后的图像可以通过对H(u,v)*F(u,v)得到（将负值当0灰度级显示）

滤波在频率域中更加直观，在空间域中滤波器模板越小越好

4.3实现图像平滑的频率滤波器

由于边缘和其他尖锐变化（如噪声）在图像灰度级中主要处于傅里叶变化的高频部分，因此平滑（模糊）可以通过衰减指定图像傅里叶变化中高频成分的范围来实现。
频率域平滑滤波器有
1、理想低通滤波器

2、巴特沃思低通滤波器

3、高斯低通滤波器

√ 边缘和噪声等尖锐变化处于傅里叶变换的高 频部分

√ 平滑可以通过衰减高频成分的范围来实现

√ 理想低通滤波器：尖锐

√ 巴特沃思低通滤波器：处于理想和高斯滤波 器之间

√ 高斯低通滤波器：平滑

1.理想低通滤波器（ILPF）

在以原点为圆心，以D0为半径的圆内，无衰减的通过所有频率，而在该圆外切断所有频率的二维低通滤波器
理想低通滤波器是最简单的低通滤波器，是通过截断傅里叶变换中所有的高频成分，这些高频成分距离变换原点的距离比指定距离D0要远的多。
其变换函数是：
在这里插入图片描述
在这里插入图片描述
在这里插入图片描述
理想低通滤波器有振铃和模糊现象
理想低通滤波器在频率域的形状为矩形，那么其傅立叶逆变换在时间域为sinc函数；
图像处理中，对一幅图像进行滤波处理，若选用的频域滤波器具有陡峭的变化，则会使滤波图像产生“振铃”，所谓“振铃”，就是指输出图像的灰度剧烈变化处产生的震荡，就好像钟被敲击后产生的空气震荡

实现代码：

#include<opencv2/opencv.hpp>
#include<iostream>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = {  src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数
	
}


//定义理想低通滤波器函数
Mat ideal_low_kernel(Mat&src, float sigma)
{
	Mat ideal_low_pass(src.size(), CV_32FC1);
	float d0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			double d = sqrt(pow((i - src.rows / 2), 2) + pow((j - src.cols / 2), 2));//分子,计算pow必须为float型
			if (d <= d0)
			{
				ideal_low_pass.at<float>(i, j) = 1;
			}
			else
			{
				ideal_low_pass.at<float>(i, j) = 0;
			}
		}
	}
	//显示滤波器
	string name = "理想低通滤波器=" + std::to_string(sigma);//std::to_string将数值转化为字符串。返回对应的字符串。
	cv::imshow(name, ideal_low_pass);
	return ideal_low_pass;
}
cv::Mat ideal_low_pass_filter(Mat &src, float sigma)
{
	Mat padded = image_make_border(src);
	Mat ideal_kernel = ideal_low_kernel(padded, sigma);
	Mat result = frequence_filter(padded, ideal_kernel);
	return result;
}
int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal = ideal_low_pass_filter(img, 10.00000);
	ideal = ideal(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img",ideal);
	waitKey(0);
	return 0;
}

D0=10时结果：
在这里插入图片描述

2.巴特沃斯低通滤波器（BLPF）

截止频率位于距原点D0处的n阶巴特沃斯低通滤波器（BLPF）的传递函数为：
在这里插入图片描述
与理想滤波器不同的是，BLPF传递函数并没有给出明显截止的尖锐的不连续性。
在上式中，截止频率点是当D(u，v)=D0时的点【即H(u，v)从其最大值下降为50%】
空间域的一阶巴特沃斯滤波器没有振铃现象，在二阶滤波器中，振铃现象通常很难察觉，只有更高阶数才能显示明显的振铃现象。
//巴特沃斯低通滤波器

#include<opencv2/opencv.hpp>
#include<iostream>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = {  src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数
	
}
Mat butterworth_low_kernel(Mat &scr, float sigma, int n)//n是阶数
{
	Mat butterworth_low_pass(scr.size(), CV_32FC1); //，CV_32FC1
	double D0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for(int i = 0; i < scr.rows; i++)
	{
		for (int j = 0; j < scr.cols; i++)
		{
			//下面两行是应用巴特沃斯的公式
			double d = sqrt(pow((i - scr.rows / 2), 2) + pow((j - scr.cols / 2), 2));//分子,计算pow必须为float型
			butterworth_low_pass.at<float>(i, j) = 1.0 / (1 + pow(d/D0, 2 * n));
		}
	}
	string name = "巴特沃斯低通滤波器d0=" + std::to_string(sigma) + "n=" + std::to_string(n);
	imshow(name, butterworth_low_pass);
	return butterworth_low_pass;
}
//巴特沃斯函数封装
Mat butterworth_low_paass_filter(Mat &src, float d0, int n)
{
	//H = 1 / (1+(D/D0)^2n)    n表示巴特沃斯滤波器的次数
   //阶数n=1 无振铃和负值    阶数n=2 轻微振铃和负值  阶数n=5 明显振铃和负值   阶数n=20 与ILPF相似
	Mat padded = image_make_border(src);
	Mat butterworth_kernel = butterworth_low_kernel(padded, d0, n);
	Mat result = frequence_filter(padded, butterworth_kernel);
	return result;
}

3.高斯低通滤波器GLPF

高斯滤波器没有振铃现象，且边缘平滑，传递函数为：其中D0是截止频率，当D=D0时，GLPF下降到其最大值的0.607处。
在这里插入图片描述
在这里插入图片描述

#include<opencv2/opencv.hpp>
#include<iostream>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = {  src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数	
}
//高斯低通滤波器
Mat gaussian_low_pass_kernel(Mat scr, float sigma)
{
	Mat gaussianBlur(scr.size(), CV_32FC1); //，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型）
		}
	}
	string name = "高斯低通滤波器d0=" + std::to_string(sigma);
	imshow(name, gaussianBlur);
	return gaussianBlur;
}
Mat gaussian_low_pass_filter(Mat &src, float d0)
{
	Mat padded = image_make_border(src);
	Mat gaussian_kernel = gaussian_low_pass_kernel(padded, d0);//理想低通滤波器
	Mat result = frequence_filter(padded, gaussian_kernel);
	return result;
}
int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal = gaussian_low_pass_filter(img, 160);
	ideal = ideal(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img",ideal);
	waitKey(0);
	return 0;

}

D0=160时结果
在这里插入图片描述

4.4使用频率域滤波器锐化图像

高通滤波衰减傅里叶变换中的低频分量而不影响高频分量。
理想高通滤波器：锐利性大
巴特沃斯高通滤波器：介于上下两者之间，属于过渡
高斯高通滤波器：宽阔平滑性大

4.4.1理想高通滤波器(二维IHPF)

二维理想高通滤波器与理想低通滤波器相对，传递函数为：
在这里插入图片描述
IHPF是将以D0为半径的圆内的所有频率置零，而毫无衰减的通过圆外所有频率。和理想低通一样IHPF在物理上是不可实现的。
IHPF的振铃现象很严重，导致失真。
所有图像的恒定背景在高通滤波后的图像中都是0，因此高通滤波类似于空间域的差分。
典型理想高通滤波器的透视图、图像表示、剖面图
典型的频率域理想高通滤波器的空间表示及通过滤波器中心的对应灰度剖面图
实现代码：

#include<opencv2/opencv.hpp>
#include<iostream>
#include<string>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = { src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数

}
//理想高通滤波器传递函数的实现
Mat ideal_high_kernel(Mat &src, float sigma)
{
	Mat ideal_high_pass(src.size(), CV_32FC1);
	float d0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			double d = sqrt(pow((i - src.rows / 2), 2) + pow((j - src.cols / 2), 2));
			if (d <= d0)
			{
				ideal_high_pass.at<float>(i, j) = 0;
			}
			else
			{
				ideal_high_pass.at<float>(i, j) = 1;
			}
		}
	}
	string name = "理想高通滤波器d0=" + std::to_string(sigma);
	imshow(name, ideal_high_pass);
	return ideal_high_pass;
}
//对理想高通滤波器函数进行封装
cv::Mat ideal_high_pass_filter(Mat &src, float sigma)
{
	Mat padded = image_make_border(src);
	Mat ideal_kernel = ideal_high_kernel(padded, sigma);
	Mat result = frequence_filter(padded, ideal_kernel);
	return result;
}
int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal = ideal_high_pass_filter(img, 160);
	ideal = ideal(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img", ideal);
	waitKey(0);
	return 0;

}

D0=160的结果：
在这里插入图片描述

4.4.2巴特沃斯高通滤波器BHPF

截止频率为D0的n阶巴特沃斯高通滤波器定义为：
在这里插入图片描述
其中D(u，v)：
在这里插入图片描述
截止频率越大，使用BHPF得到的结果越平滑
在这里插入图片描述
在这里插入图片描述

#include<opencv2/opencv.hpp>
#include<iostream>
#include<string>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = { src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数

}
//巴特沃斯高通滤波器
cv::Mat butterworth_high_kernel(Mat &src, float sigma, int n)
{
	Mat butterworth_high_pass(src.size(), CV_32FC1);
	double d0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			double d = sqrt(pow((i - src.rows / 2), 2) + pow((j - src.cols / 2), 2));
			butterworth_high_pass.at<float>(i, j) = 1.0 / (1 + pow(d0 / d, 2 * n));
		}
	}
	string name = "巴特沃斯高通滤波器d0=" + std::to_string(sigma) + "n = " + std::to_string(n);
	imshow(name, butterworth_high_pass);
	return butterworth_high_pass;
}
//巴特沃斯高通滤波器封装
Mat butterworth_high_pass_filter(Mat &src, float d0, int n)
{
	// H = 1 / (1 + (D0 / D) ^ 2n)    n表示巴特沃斯滤波器的次数
	Mat padded = image_make_border(src);
	Mat butterworth_kernel = butterworth_high_kernel(padded, d0, n);
	Mat result = frequence_filter(padded, butterworth_kernel);
	return result;
}
int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal2 = butterworth_high_pass_filter(img, 160,2);
	ideal2 = ideal2(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img", ideal2);
	waitKey(0);
	return 0;

}

结果显示：
在这里插入图片描述

4.4.3高斯高通滤波器GHPF

对微小物体和细线条采用高斯滤波，结果也比较清晰
截止频率处在距频率矩形中心距离为D0的高斯高通滤波器传递函数：
在这里插入图片描述
在这里插入图片描述

#include<opencv2/opencv.hpp>
#include<iostream>
#include<string>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = { src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数
}
//高斯高通滤波器传递函数
Mat gaussian_high_pass_kernel(Mat &src, float sigma)
{
	Mat gaussianBlur(src.size(), CV_32FC1);
	float d0 = 2 * sigma*sigma;//高斯函数参数，越小，频率高斯滤波器越窄，滤除高频成分越多，图像就越平滑
	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			float d = sqrt(pow(float(i - src.rows / 2), 2) + pow(float(j - src.cols / 2), 2));
			gaussianBlur.at<float>(i, j) = 1 - expf(-pow(d,2) /2*pow(d0,2));

		}
	}
	string name = "高斯高通滤波器d0=" + std::to_string(sigma);
	imshow(name,gaussianBlur);
	return gaussianBlur;
}
//高斯高通滤波器封装
Mat gaussian_high_pass_filter(Mat &src, float d0)
{
	Mat padded = image_make_border(src);
	Mat gaussian_kernel = gaussian_high_pass_kernel(padded, d0);//理想低通滤波器
	Mat result = frequence_filter(padded, gaussian_kernel);
	return result;

}

int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal3 = gaussian_high_pass_filter(img, 160);
	ideal3 = ideal3(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img3", ideal3);
	waitKey(0);
	return 0;
}

在这里插入图片描述

4.4.4频率域中的拉普拉斯算子

首先是拉普拉斯在频率域中使用的滤波器模板如下：
在这里插入图片描述
或者采用下面在这里插入图片描述
![![](https://img-blog.csdnimg.cn/2021052011031447.png)

然后计算拉普拉斯算子图像（F(U,V）是f的傅里叶变换：
在这里插入图片描述
最后与源图像叠加，实现源图像的增强：

在这里插入图片描述
在频率域中上述式子可写成：
在这里插入图片描述

#include<opencv2/opencv.hpp>
#include<iostream>
#include<string>

using namespace std;
using namespace cv;

//首先对图片定义一个函数，获取傅里叶最优尺寸大小，进而扩充图像边界，避免缠绕
cv::Mat image_make_border(cv::Mat &src)
{
	int w = getOptimalDFTSize(src.cols);//getOptimalDFTSize函数返回给定向量尺寸的傅里叶最优尺寸大小
	int h = getOptimalDFTSize(src.rows);
	Mat padded;
	copyMakeBorder(src, padded, 0, h - src.rows, 0, w - src.cols, BORDER_CONSTANT, Scalar::all(0));//该函数为边缘填充函数
	padded.convertTo(padded, CV_32FC1);
	return padded;
}
//第二定义一个频率域卷积的基本公式
Mat frequence_filter(Mat& src, Mat& blur)
{
	//先进行DFT傅里叶变换
	Mat plane[] = { src, Mat::zeros(src.size() , CV_32FC1) };//创建通道，存储DFT后的实部与虚部（CV_32F，必须为单通道数）
	Mat complexIm;
	merge(plane, 2, complexIm);//合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	dft(complexIm, complexIm);//进行傅立叶变换，结果保存在自身
	//进行中心化操作
	split(complexIm, plane);//分离通道（数组分离）将原图数组分为0和1，黑的白的，分别对应实部和虚部
//    plane[0] = plane[0](Rect(0, 0, plane[0].cols & -2, plane[0].rows & -2));//这里为什么&上-2具体查看opencv文档
//    //其实是为了把行和列变成偶数 -2的二进制是11111111.......10 最后一位是0
	int cx = plane[0].cols / 2;
	int cy = plane[0].rows / 2;
	Mat part1_r(plane[0], Rect(0, 0, cx, cy)); //元素坐标表示为(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)); //元素坐标(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);

	//第三步是滤波器函数和DFT结果相乘
	Mat 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);//实部与虚部合并


	//最后得到原图频谱图
	magnitude(plane[0], plane[1], plane[0]);//获取幅度图像，0通道为实部通道，1为虚部，因为二维傅立叶变换结果是复数
	plane[0] += Scalar::all(1);//傅立叶变换后的图片不好分析，进行对数处理，结果比较好看
	log(plane[0], plane[0]);// float型的灰度空间为[0，1])
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX);//归一化便于显示

	idft(BLUR, BLUR);//idft二维傅里叶变化的结果也为复数
	split(BLUR, plane);//分离通道，主要获取通道
	magnitude(plane[0], plane[1], plane[0]); //求幅值(模)
	normalize(plane[0], plane[0], 1, 0, CV_MINMAX); //归一化便于显示
	return plane[0];//返回参数

}
//拉普拉斯滤波器传递函数
Mat laplace_kernel(Mat &src)
{
	Mat laplace_pass(src.size(), CV_32FC2);
	int row_num = src.rows;
	int col_num = src.cols;
	for (int i = 0; i < row_num; i++)
	{
		float *p = laplace_pass.ptr<float>(i);
		for (int j = 0; j < col_num; j++)
		{
			float d = pow((i - row_num / 2), 2) + pow((j - col_num / 2), 2);
			p[2 * j] = -4 * CV_PI * CV_PI * d;//应用公式
			p[2 * j + 1] = -4 * CV_PI * CV_PI * d;
		}
	}
	cv::Mat temp[] = { cv::Mat::zeros(src.size(), CV_32FC1),
					 cv::Mat::zeros(src.size(), CV_32FC1) };
	cv::split(laplace_pass, temp);
	std::string name = "laplace滤波器";
	cv::Mat show;
	cv::normalize(temp[0], show, 1, 0, CV_MINMAX);
	cv::imshow(name, show);
	return laplace_pass;
}
//laplace滤波器封装
cv::Mat laplace_filter(cv::Mat &src)
{
	cv::Mat padded = image_make_border(src);
	cv::Mat lap_kernel = laplace_kernel(padded);
	cv::Mat result =frequence_filter(padded, lap_kernel);
	return result;
}


int main()
{
	Mat img = imread("C:/Users/征途/Desktop/vs-cpp/第四章/01.jpg", IMREAD_GRAYSCALE);
	cv::Mat ideal4 = laplace_filter(img);
	ideal4 = ideal4(cv::Rect(0, 0, img.cols, img.rows));
	imshow("img4", ideal4);
	waitKey(0);
	return 0;

}

在这里插入图片描述

4.4.5 同态滤波器

同态滤波是把频率过滤和灰度变换结合起来的一种图像处理方法，它依靠图像的照度/ 反射率模型作为频域处理的基础，利用压缩亮度范围和增强对比度来改善图像的质量。使用这种方法可以使图像处理符合人眼对于亮度响应的非线性特性，避免了直接对图像进行傅立叶变换处理的失真。
1、对图片的照射分量和反射分量的处理
一幅图像f(x,y)可以表示为其照射分量i(x,y)和反射分量r(x,y)的乘积，即：
f(x,y) = i(x,y)r(x,y)
其中f(x, y)代表原图像，i(x,y)代表照射强度，r(x,y)代表反射强度，
一般的图像光照是均匀变化的，所以i应该是低频分量，而不同的物体对光的反射是具有突变型的，所以r(x,y)是高频分量
同态滤波实现步骤：
在这里插入图片描述

    上式不能直接用于对照射和反射的频率分量进行操作，因为乘积的傅里叶变换不是变换的乘积。

在这里插入图片描述
因此采用对数处理，并且针对对数式进行傅里叶变换：

在这里插入图片描述
在这里插入图片描述
或者在这里插入图片描述
在这里插入图片描述
使用一个同态滤波器H(u,v)对图像的Z(u,v)进行滤波：
在这里插入图片描述
在空间域中滤波后的图像为:
在这里插入图片描述
最后取结果的指数来消除前面所取的指数：
在这里插入图片描述
上图上面的i0(x,y)和r0(x,y)是输出（处理后）图像的照射和反射分量。
上述的方法是以称为同台系统的一类系统的特殊情况。重点是将一幅图像的照射分量和反射分量进行分离，然后再通过同态滤波H(u,v)对分量进行操作。

2、高斯同态滤波器传递函数
图像的照射分量通常由慢的空间变化来表征，而反射分量引起突变，特别是在不同物体的连接部分。
因此导致图像取对数后的傅里叶变换的低频成分与照射相联系，而高频成分与反射相联系。同态滤波器通过滤波器函数H（u,v）可以用不同的可控方法影响傅里叶变换的低频和高频分量。

高斯同态滤波器传递函数：在这里插入图片描述
其中D(U,V)是距离中心的距离，由前面的式子确定，常数c控制函数坡度的锐利度，它在γL和γH之间过度。
如果γL和γH选定， 而γL<1 和γH>1,那么下图滤波器函数趋向衰减低频（照射）的贡献，而增强高频（反射）的贡献。最终结果是同时进行动态范围的压缩和对比度的增强。
同态滤波器的传递函数一般在低频部分小于1，高频部分大于1。
在这里插入图片描述
3、实现代码如下：

#include<opencv2/opencv.hpp>
#include<iostream>
#include<string>
#include<time.h>
#include <math.h>
using namespace std;
using namespace cv;

//确定函数
Mat HomoFilter(cv::Mat src)
{
	src.convertTo(src, CV_64FC1);//转换数据类型
	int rows = src.rows;
	int cols = src.cols;
	int m = rows % 2 == 1 ? rows + 1 : rows;//三目运算符
	int n = cols % 2 == 1 ? cols + 1 : cols;
	copyMakeBorder(src, src, 0, m - rows, 0, n - cols, BORDER_CONSTANT, Scalar::all(0));//CopyMakeBorder复制图像并且制作边界。
	rows = src.rows;
	cols = src.cols;
	Mat dst(rows, cols, CV_64FC1);//输出

	//1、求对数
	for (int i = 0; i < rows; i++) {
		double *srcdata = src.ptr<double>(i);
		double *logdata = src.ptr<double>(i);
		for (int j = 0; j < cols; j++) {
			logdata[j] = log(srcdata[j] + 0.0001);
		}
	}
	//2.dct离散余弦变换
	Mat mat_dct = Mat::zeros(rows, cols, CV_64FC1);
	dct(src, mat_dct);

	//3. 高斯同态滤波器
	Mat H_u_v;
	double gammaH = 1.5;//1   高频增益
	double gammaL = 0.5;//<1   低频增益
	double C = 1;//斜面锐化常数 斜率
	double  d0 = (src.rows / 2) * (src.rows / 2) + (src.cols / 2) * (src.cols / 2);//double  d0 = 150;//5-200  截止频率 越大越亮
	double  d2 = 0;
	H_u_v = Mat::zeros(rows, cols, CV_64FC1);
	for (int i = 0; i < rows; i++)
	{
		double * dataH_u_v = H_u_v.ptr<double>(i);
		for (int j = 0; j < cols; j++) {
			d2 = pow(i, 2.0) + pow(j, 2.0);
			dataH_u_v[j] = (gammaH - gammaL) * (1 - exp(-C * d2 / d0)) + gammaL;
		}
	}
	H_u_v.ptr<double>(0)[0] = 1.1;
	mat_dct = mat_dct.mul(H_u_v);

	//4.idct
	idct(mat_dct, dst);//离散余弦反转换的公式
	//exp
	for (int i = 0; i < rows; i++) 
	{
		double  *srcdata = dst.ptr<double>(i);
		double *dstdata = dst.ptr<double>(i);
		for (int j = 0; j < cols; j++)
		{
			dstdata[j] = exp(srcdata[j]);
		}
	}
	dst.convertTo(dst, CV_8UC1);
	return dst;

}

//主函数
int main()
{
	Mat src = cv::imread("C:/Users/征途/Desktop/vs-cpp/第四章/同态.jpg");
	imshow("origin", src);
	int originrows = src.rows;
	int origincols = src.cols;
	Mat dst(src.rows, src.cols, CV_8UC3);
	cvtColor(src, src, COLOR_BGR2YUV);//转化为YUV通道

	vector <Mat> yuv;//创建vector对象yuv
	split(src, yuv);//将src通道分离
	Mat nowY = yuv[0];
	Mat newY = HomoFilter(nowY);
	Mat tempY(originrows, origincols, CV_8UC1);
	for (int i = 0; i < originrows; i++) 
	{
		for (int j = 0; j < origincols; j++) 
		{
			tempY.at<uchar>(i, j) = newY.at<uchar>(i, j);
		}
	}
	yuv[0] = tempY;
	merge(yuv, dst);//通道合并，merge()是C++标准库的函数，主要实现函数的排序和合并
	cvtColor(dst, dst, COLOR_YUV2BGR);
	imshow("result", dst);
	waitKey(0);
}

在这里插入图片描述

cv吴彦祖
关注
  • 9
    点赞
  • 踩
  • 85
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
图像滤波----低通滤波，中值滤波，高通滤波，方向滤波（Sobel），拉普拉斯变换
渣渣的博客
11-24 3万+
①观察灰度分布来描述一幅图像成为空间域，观察图像变化的频率被成为频域。 ②频域分析：低频对应区域的图像强度变化缓慢，高频对应的变化快。低通滤波器去除了图像的高频部分，高通滤波器去除了图像的低频部分。(1)低通滤波 ①栗子：#include <iostream> #include <opencv2/core/core.hpp> #include <opencv2/imgproc/imgproc.h
利用低通滤波器和陷波滤波器对数字图像进行频域滤波
green_building_xliu的博客
04-14 7459
1. 对图像做傅里叶变换，得到其频域图像 2 设计理想低通滤波器 3.滤波并显示结果图像  4.设计陷波滤波器 5 滤波显示结果图像查看不同参数的效 matlab代码： clear all; close all; clc; %1.对图像做傅里叶变换，得到其频域图像 src = im2double(imread('image3.bmp')); src = rgb2g
频域滤波器
nameofcsdn的博客
06-11 2127
。。。
无处不在的拉普拉斯——边缘，斑点，金字塔
zcg1942的博客
08-23 930
作为拿破仑的老师，拉普拉斯有 拉普拉斯算子 二阶导的过零点可以用来检测边缘 有限差分算子，可以由二阶导定义或者泰勒展开推导出。 高斯拉普拉斯LoG 因为二阶导对噪声敏感，所以通过先使用高斯模糊，而我们可以提前计算高斯滤波的拉普拉斯响应。画出图像是墨西哥草帽。 墨西哥草帽和图像做相关就可以实现斑点检测（结合极大值的寻找） LoG可以使用DoG近似。DoG是不同的尺度因子，所以两个高斯滤波核中心会错位，二者做差就类似拉普拉斯中的左右或者上下做差。 拉普拉斯金字塔 先下采样再上采样就会损失信息，
高斯函数及高斯滤波器
热门推荐
贾路飞的博客
07-02 3万+
1 一维高斯分布 1.1 一维高斯分布的定义 若连续型随机变量X的概率密度为： 其中，为常数，则称X服从参数为，的正态分布或高斯分布，记为 1.2 一维高斯分布的曲线 横轴表示可能的取值x，竖轴表示概率分布密度F(x)，那么不难理解这样一个曲线与x轴围...
《OpenCv视觉之眼》Python图像处理七 :Opencv图像处理之高通滤波和低通滤波原理及构造
qq_42451251的博客
08-07 6539
本专栏主要介绍如果通过OpenCv-Python进行图像处理，通过原理理解OpenCv-Python的函数处理原型，在具体情况中，针对不同的图像进行不同等级的、不同方法的处理，以达到对图像进行去噪、锐化等一系列的操作。同时，希望观看本专栏的小伙伴可以理解到OpenCv进行图像处理的强大哦，如有转载，请注明出处(原文链接和作者署名)，感谢各位小伙伴啦！ 前文参考： 《OpenCv视觉之眼》Python图像处理一 :Opencv-python的简介及Python环境搭建 《OpenCv视觉之眼》Python图像
MATLAB中滤波函数、频率响应函数以及频率响应函数不同表达形式的转换
yindq1220的博客
07-11 8588
MATLAB中滤波函数filter、频率响应函数freqz以及频率响应函数不同表达形式的转换tf2zp/zp2tf
使用同态滤波+拉普拉斯算法处理图像(APMCM为例)
12-22
同态滤波 该方法在博客中有详细说明。主要思想是： 图像f(x,y)f(x,y)f(x,y)可由fi(x,y)fr(x,y)f_i(x,y)f_r(x,y)fi​(x,y)fr​(x,y)组成，然后两边取对数并做傅里叶变换后就有 DFT[lnf(x,y)]=DFT[lnfi(x,y)]+DFT[lnfr...
matlab图像技术：10 高通滤波、高斯滤波、拉普拉斯算子锐化、理想高通滤波、同态滤波
06-28
本压缩包文件“matlab图像技术：10 BUTTERWORTH高通滤波、高斯滤波、拉普拉斯算子锐化、理想高通滤波、同态滤波、伪彩色增强实现图像增强.zip”提供了多个图像处理技术的实现，让我们逐一探讨这些关键技术。...
数字图像处理第三版（Gonzalez）复习整理
03-09
《数字图像处理第三版（Gonzalez）复习整理》由Jacques Sun师兄精心整理，涵盖了图像处理领域的核心概念，包括像素邻域、邻接性、距离度量、灰度变换、直方图处理、滤波器以及傅里叶变换等多个方面。以下是对这些...
图像处理中的低通滤波器
03-09
程序中为图像处理的各种低通滤波器，主要包括理想低通滤波器、梯形低通滤波器、巴特沃斯低通滤波器、指数低通滤波器.
巴特沃斯低通滤波器处理彩色图像代码
12-17
可以运行，本人的大作业，处理的彩色图片，对三通道分别处理
matlab图像专题：41 同态滤波实现图像增强.zip
最新发布
08-23
在同态滤波中，我们通常会使用一个低通滤波器来处理线性部分，以保留图像的基本结构；同时，使用一个高通滤波器来处理非线性部分，以去除噪声。 接下来，我们将讨论如何构建同态滤波的步骤： 1. **预处理**：加载...
《数字图像处理》冈萨雷斯 中文+英文 pdf
12-17
全书分为12章，包括绪论，数字图像基础，灰度变换与空间滤波，频域滤波，图像复原与重建，彩色图像处理，小波及多分辨率处理，图像压缩，形态学图像处理，图像分割，表现与描述，目标识别。 目录编辑 前言15 致谢19 ...
自适应滤波器 | 频域自适应滤波器（FLMS）
qq_45490227的博客
10-11 2713
频域自适应滤波器算法实现
灰度图像的频率域滤波——拉普拉斯高通滤波（Matlab）
lengo的专栏
09-03 2097
clc; clear all; close all; %读入图像，并转换为double型 I=imread('D:\Gray Files\4-41.tif'); I_D=im2double(I); %获得图像的高度和宽度 [M,N]=size(I_D); %图像中心点 M0=M/2; N0=N/2; J=fft2(I_D); J_shift=fftshift(J); %%%=========...
图卷积（6）——过平滑现象（1）
maple的博客
11-07 5346
图卷积（6）——过平滑现象（1） 过平滑现象的提出是在论文Deeper Insights into Graph Convolutional Networks for Semi-Supervised Learning中。作者将GCN用于半监督学习任务，发现采用多层GCN之后，会出现过平滑现象（over-smoothing）。 文章提出三个核心观点： 拉普拉斯平滑：GCN中图卷积是特殊形式的拉普拉斯平滑，可以聚合临近节点的特征，使同一类别的节点的特征相似； 过平滑问题：多层GCN使得节点输出特征过度平滑，这样
OpenCV：条纹图像频谱分离
weixin_43865692的博客
08-12 1798
条纹函数及条纹傅里叶变换 在生成条纹图像时，用到余弦函数为：m+m*cos(2*π*f) 而其傅里叶变换为： 具体条纹图像生成代码参考：https://blog.csdn.net/weixin_43865692/article/details/119141967 从傅里叶变换的函数可以看出其条纹频谱有三个频率，分别为中心点以及关于中心对称的两点。 条纹图像频谱分离原理 在此选用三个等式求得三个频率。选取三个不同相移的余弦函数，相移分别为0，2*π/3，4*π/3。再得到对应的傅里叶变换的函数
傅里叶变换、拉普拉斯变换、z变换之间的联系
weixin_45024950的博客
09-05 5452
本文阐述了傅里叶变换、拉普拉斯变换、z变换之间的联系
图像低通滤波和高通滤波
06-07
图像滤波是一种常见的图像处理方法，可以用来去除图像中的噪声、增强图像的边缘等。图像滤波可以分为低通滤波和高通滤波两种。 低通滤波可以去除图像中高频部分，保留图像中低频部分，实现图像平滑处理。低通滤波器通常被用来去除图像中的高频噪声，使得图像更加清晰，常见的低通滤波器有均值滤波、高斯滤波等。 高通滤波可以增强图像中高频部分，抑制图像中低频部分，实现图像锐化处理。高通滤波器通常被用来增强图像的边缘和细节，常见的高通滤波器有拉普拉斯滤波、Sobel滤波、Prewitt滤波等。
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值