OpenCV-一维频域滤波器

03669748

作者：翟天保Steven
版权声明：著作权归作者所有，商业转载请联系作者获得授权，非商业转载请注明出处

场景需求

       之前分享了许多图像频域滤波器的功能，如理想滤波器、高斯滤波器、巴特沃斯滤波器等等，都是二维形式的写法；有粉丝提出，不知道一维数据该怎么进行类似的处理。其实搞懂原理，二维降一维也是非常方便的，为了方便大家使用，本文将用C++和OpenCV实现一维频域滤波器功能。

       示例为巴特沃斯，如果想换成别的，只要把滤波核更换即可。

相关功能函数的C++实现代码

       二维与一维的不同之处：

  1. 获取DFT变换最佳高度时，一维只需要考虑列或者行，二维需要都考虑。
  2. 频谱搬移操作，二维是四象限，一维是两部分，即上下或者左右交换即可。
  3. 频域中滤波器尺寸与距离中心的距离有关，二维的距离公式是一个圆表达式，一维的距离公式是差值绝对值。
// 图像边界处理（一维）
cv::Mat image_make_border_onedim(cv::Mat &src)
{
	int h = cv::getOptimalDFTSize(src.rows); // 获取DFT变换的最佳高度

	cv::Mat padded;
	// 常量法扩充图像边界，常量 = 0
	cv::copyMakeBorder(src, padded, 0, h - src.rows, 0, 0, cv::BORDER_CONSTANT, cv::Scalar::all(0));
	padded.convertTo(padded, CV_32FC1);

	return padded;
}

// fft变换后进行频谱搬移（一维）
void fftshift_onedim(cv::Mat &plane0, cv::Mat &plane1)
{
	// 以下的操作是移动图像  (零频移到中心)
	int cy = plane0.rows / 2;
	cv::Mat part1_r(plane0, cv::Rect(0, 0, 1, cy));  // 元素坐标表示为(cx, cy)
	cv::Mat part2_r(plane0, cv::Rect(0, cy, 1, cy));

	cv::Mat temp;
	part1_r.copyTo(temp);  //上与下交换位置(实部)
	part2_r.copyTo(part1_r);
	temp.copyTo(part2_r);

	cv::Mat part1_i(plane1, cv::Rect(0, 0, 1, cy));  //元素坐标(cx,cy)
	cv::Mat part2_i(plane1, cv::Rect(0, cy, 1, cy));

	part1_i.copyTo(temp);  //上与下交换位置(虚部)
	part2_i.copyTo(part1_i);
	temp.copyTo(part2_i);
}

// pow操作
Mat powZ(cv::InputArray src, double power) {
	cv::Mat dst;
	cv::pow(src, power, dst);
	return dst;
}

// sqrt操作
Mat sqrtZ(cv::InputArray src) {
	cv::Mat dst;
	cv::sqrt(src, dst);
	return dst;
}

// 频率域滤波
cv::Mat frequency_filter(cv::Mat &scr, cv::Mat &blur)
{
	cv::Mat mask = scr == scr;
	scr.setTo(0.0f, ~mask);

	//创建通道，存储dft后的实部与虚部（CV_32F，必须为单通道数）
	cv::Mat plane[] = { scr.clone(), cv::Mat::zeros(scr.size() , CV_32FC1) };

	cv::Mat complexIm;
	cv::merge(plane, 2, complexIm); // 合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	cv::dft(complexIm, complexIm); // 进行傅立叶变换，结果保存在自身

	// 分离通道（数组分离）
	cv::split(complexIm, plane);

	// 以下的操作是频域迁移
	fftshift_onedim(plane[0], plane[1]);

	// *****************滤波器函数与DFT结果的乘积****************
	cv::Mat blur_r, blur_i, BLUR;
	cv::multiply(plane[0], blur, blur_r);  // 滤波（实部与滤波器模板对应元素相乘）
	cv::multiply(plane[1], blur, blur_i);  // 滤波（虚部与滤波器模板对应元素相乘）
	cv::Mat plane1[] = { blur_r, blur_i };

	// 再次搬移回来进行逆变换
	fftshift_onedim(plane1[0], plane1[1]);
	cv::merge(plane1, 2, BLUR); // 实部与虚部合并

	cv::idft(BLUR, BLUR);       // idft结果也为复数
	BLUR = BLUR / BLUR.rows / BLUR.cols;

	cv::split(BLUR, plane);//分离通道，主要获取通道

	return plane[0];
}

// 巴特沃斯低通滤波核函数（一维）
cv::Mat butterworth_low_kernel_onedim(cv::Mat &scr, float sigma, int n)
{
	cv::Mat butterworth_low_pass(scr.size(), CV_32FC1); //，CV_32FC1
	float D0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for (int i = 0; i < scr.rows; i++) {
		float d = abs(float(i - scr.rows / 2));
		butterworth_low_pass.at<float>(i, 0) = 1.0f / (1.0f + pow(d / D0, 2 * n));
	}
	return butterworth_low_pass;
}

// 巴特沃斯低通滤波（一维）
cv::Mat butterworth_low_pass_filter_onedim(cv::Mat &src, float d0, int n)
{
	// H = 1 / (1+(D/D0)^2n)   n表示巴特沃斯滤波器的次数
	// 阶数n=1 无振铃和负值    阶数n=2 轻微振铃和负值  阶数n=5 明显振铃和负值   阶数n=20 与ILPF相似
	cv::Mat padded = image_make_border_onedim(src);
	cv::Mat butterworth_kernel = butterworth_low_kernel_onedim(padded, d0, n);
	cv::Mat result = frequency_filter(padded, butterworth_kernel);
	return result;
}

测试代码

#include <iostream>
#include <fstream>
#include <sstream>
#include <vector>
#include<opencv2/opencv.hpp>
#include<ctime>

using namespace cv;
using namespace std;

// 读取数据
cv::Mat ReadOneDimData(string &file)
{
	cv::Mat result;
	std::ifstream in(file);
	std::string str;

	getline(in, str);
	while (str != "")
	{
		int a = str.find(',');
		string front = str.substr(0, a);
		string back = str.substr(a + 1);
		float wave = stof(front);
		float h = stof(back);
		cv::Mat temp = (cv::Mat_<float>(1, 2) << wave, h);
		result.push_back(temp);
		getline(in, str);
	}

	return result;
}

// 显示数据曲线
cv::Mat ShowDataCurve(cv::Mat input)
{
	// 归一化
	cv::Mat curve = input.clone();
	cv::normalize(curve, curve, 0, 255, cv::NORM_MINMAX);
	curve.convertTo(curve, CV_8UC1);

	// 绘制简易曲线图
	int size = input.rows;
	cv::Mat pic = cv::Mat::zeros(256, size, CV_8UC1);
	for (int i = 0; i < size; ++i)
	{
		uchar h = curve.at<uchar>(i, 0);
		pic.at<uchar>(255 - h, i) = 255;
	}
	return pic;
}

// 图像边界处理（一维）
cv::Mat image_make_border_onedim(cv::Mat &src)
{
	int h = cv::getOptimalDFTSize(src.rows); // 获取DFT变换的最佳高度

	cv::Mat padded;
	// 常量法扩充图像边界，常量 = 0
	cv::copyMakeBorder(src, padded, 0, h - src.rows, 0, 0, cv::BORDER_CONSTANT, cv::Scalar::all(0));
	padded.convertTo(padded, CV_32FC1);

	return padded;
}

// fft变换后进行频谱搬移（一维）
void fftshift_onedim(cv::Mat &plane0, cv::Mat &plane1)
{
	// 以下的操作是移动图像  (零频移到中心)
	int cy = plane0.rows / 2;
	cv::Mat part1_r(plane0, cv::Rect(0, 0, 1, cy));  // 元素坐标表示为(cx, cy)
	cv::Mat part2_r(plane0, cv::Rect(0, cy, 1, cy));

	cv::Mat temp;
	part1_r.copyTo(temp);  //上与下交换位置(实部)
	part2_r.copyTo(part1_r);
	temp.copyTo(part2_r);

	cv::Mat part1_i(plane1, cv::Rect(0, 0, 1, cy));  //元素坐标(cx,cy)
	cv::Mat part2_i(plane1, cv::Rect(0, cy, 1, cy));

	part1_i.copyTo(temp);  //上与下交换位置(虚部)
	part2_i.copyTo(part1_i);
	temp.copyTo(part2_i);
}

// pow操作
Mat powZ(cv::InputArray src, double power) {
	cv::Mat dst;
	cv::pow(src, power, dst);
	return dst;
}

// sqrt操作
Mat sqrtZ(cv::InputArray src) {
	cv::Mat dst;
	cv::sqrt(src, dst);
	return dst;
}

// 频率域滤波
cv::Mat frequency_filter(cv::Mat &scr, cv::Mat &blur)
{
	cv::Mat mask = scr == scr;
	scr.setTo(0.0f, ~mask);

	//创建通道，存储dft后的实部与虚部（CV_32F，必须为单通道数）
	cv::Mat plane[] = { scr.clone(), cv::Mat::zeros(scr.size() , CV_32FC1) };

	cv::Mat complexIm;
	cv::merge(plane, 2, complexIm); // 合并通道 （把两个矩阵合并为一个2通道的Mat类容器）
	cv::dft(complexIm, complexIm); // 进行傅立叶变换，结果保存在自身

	// 分离通道（数组分离）
	cv::split(complexIm, plane);

	// 以下的操作是频域迁移
	fftshift_onedim(plane[0], plane[1]);

	// *****************滤波器函数与DFT结果的乘积****************
	cv::Mat blur_r, blur_i, BLUR;
	cv::multiply(plane[0], blur, blur_r);  // 滤波（实部与滤波器模板对应元素相乘）
	cv::multiply(plane[1], blur, blur_i);  // 滤波（虚部与滤波器模板对应元素相乘）
	cv::Mat plane1[] = { blur_r, blur_i };

	// 再次搬移回来进行逆变换
	fftshift_onedim(plane1[0], plane1[1]);
	cv::merge(plane1, 2, BLUR); // 实部与虚部合并

	cv::idft(BLUR, BLUR);       // idft结果也为复数
	BLUR = BLUR / BLUR.rows / BLUR.cols;

	cv::split(BLUR, plane);//分离通道，主要获取通道

	return plane[0];
}

// 巴特沃斯低通滤波核函数（一维）
cv::Mat butterworth_low_kernel_onedim(cv::Mat &scr, float sigma, int n)
{
	cv::Mat butterworth_low_pass(scr.size(), CV_32FC1); //，CV_32FC1
	float D0 = sigma;//半径D0越小，模糊越大；半径D0越大，模糊越小
	for (int i = 0; i < scr.rows; i++) {
		float d = abs(float(i - scr.rows / 2));
		butterworth_low_pass.at<float>(i, 0) = 1.0f / (1.0f + pow(d / D0, 2 * n));
	}
	return butterworth_low_pass;
}

// 巴特沃斯低通滤波（一维）
cv::Mat butterworth_low_pass_filter_onedim(cv::Mat &src, float d0, int n)
{
	// H = 1 / (1+(D/D0)^2n)   n表示巴特沃斯滤波器的次数
	// 阶数n=1 无振铃和负值    阶数n=2 轻微振铃和负值  阶数n=5 明显振铃和负值   阶数n=20 与ILPF相似
	cv::Mat padded = image_make_border_onedim(src);
	cv::Mat butterworth_kernel = butterworth_low_kernel_onedim(padded, d0, n);
	cv::Mat result = frequency_filter(padded, butterworth_kernel);
	return result;
}

int main(void)
{
	// 读取数据
	string datafile = "data4.txt";
	cv::Mat data = ReadOneDimData(datafile);

	// 获取原始曲线
	cv::Mat test = data.col(1).clone();
	cv::Mat show1 = ShowDataCurve(test);

	// 获取低通数据曲线
	float D0 = 50.0f;
	Mat lowpass = butterworth_low_pass_filter_onedim(test, D0, 2);
	cv::Mat show2 = ShowDataCurve(lowpass);

	imshow("curve1", show1);
	imshow("curve2", show2);
	waitKey(0);

	system("pause");
	return 0;
}

测试效果

图1 原始数据
图2 滤波处理后数据

       为了方便大家看出差异，我写了一个简单的数据曲线显示函数，将就着看哈~上述示例数据为光谱仪的一组数据，原始数据是高频变化的信息，经过低通滤波器处理后，保留低频成分，过滤掉高频成分，因此曲线变得平滑了。

       另外，如果我的代码有什么问题，欢迎大家提出异议批评指正，一同进步~

       如果文章帮助到你了，可以点个赞让我知道，我会很快乐~加油！

© 版权声明
文章版权归作者所有，未经允许请勿转载。
THE END
综合
# OpenCV详细教程
喜欢就支持一下吧
点赞748 分享
翟天保Steven的头像-卡核
OpenCV-向上采样cv::pyrUp&&向下采样cv::pyrDown-卡核
OpenCV-向上采样cv::pyrUp&&向下采样cv::pyrDown
OpenCV-向上采样cv::pyrUp&&向下采样cv::pyrDown
2年前 1W+
OpenCV-基于自定义色条实现灰度图上色-卡核
OpenCV-基于自定义色条实现灰度图上色
OpenCV-基于自定义色条实现灰度图上色
2年前 1W+
OpenCV-文字绘制cv::putText-卡核
OpenCV-文字绘制cv::putText
OpenCV-文字绘制cv::putText
2年前 9993
OpenCV-Laplacian边缘检测-卡核
OpenCV-Laplacian边缘检测
OpenCV-Laplacian边缘检测
2年前 9971
OpenCV-巴特沃斯低通&高通滤波器（C++）-卡核
OpenCV-巴特沃斯低通&高通滤波器（C++）
OpenCV-巴特沃斯低通&高通滤波器（C++）
2年前 9957
OpenCV-寻找轮廓cv::findContours&绘制轮廓cv::drawContours-卡核
OpenCV-寻找轮廓cv::findContours&绘制轮廓cv::drawContours
OpenCV-寻找轮廓cv::findContours&绘制轮廓cv::drawContours
2年前 9840
相关推荐
排查软件关闭时访问了0xfeeefeee内存地址导致内存访问违例的崩溃-卡核
排查软件关闭时访问了0xfeeefeee内存地址导致内存访问违例的崩溃
排查软件关闭时访问了0xfeeefeee内存地址导致内存访问违例的崩溃
2个月前 4.4W+
引发0xC0000005内存违例几种可能原因分析-卡核
引发0xC0000005内存违例几种可能原因分析
引发0xC0000005内存违例几种可能原因分析
2个月前 4.3W+
卡核知识库-卡核
卡核知识库
卡核知识库
9个月前 4.2W+
使用chatGPT,我做了一些好玩的功能，内附免费使用chatGPT网址。-卡核
使用chatGPT,我做了一些好玩的功能，内附免费使用chatGPT网址。
使用chatGPT,我做了一些好玩的功能，内附免费使用chatGPT网址。
9个月前 4W+
【侯捷】C++标准11-14-卡核
【侯捷】C++标准11-14
【侯捷】C++标准11-14
4天前 3.9W+
Radioss干货丨一次冲击不够爽，连续加载考虑一下吗？-卡核
Radioss干货丨一次冲击不够爽，连续加载考虑一下吗？
Radioss干货丨一次冲击不够爽，连续加载考虑一下吗？
2年前 3.8W+
评论 抢沙发
头像
欢迎您留下宝贵的见解！
提交
头像

昵称

取消
昵称

请填写用户信息：

  • 昵称(必填)
  • 邮箱(必填)
表情
[aoman][baiyan][bishi][bizui][cahan][ciya][dabing][daku][deyi][doge][fadai][fanu][fendou][ganga][guzhang][haixiu][hanxiao][zuohengheng][zhuakuang][zhouma][zhemo][zhayanjian][zaijian][yun][youhengheng][yiwen][yinxian][xu][xieyanxiao][xiaoku][xiaojiujie][xia][wunai][wozuimei][weixiao][weiqu][tuosai][tu][touxiao][tiaopi][shui][se][saorao][qiudale][qinqin][qiaoda][piezui][penxue][nanguo][liulei][liuhan][lenghan][leiben][kun][kuaikule][ku][koubi][kelian][keai][jingya][jingxi][jingkong][jie][huaixiao][haqian][aini][OK][qiang][quantou][shengli][woshou][gouyin][baoquan][aixin][bangbangtang][xiaoyanger][xigua][hexie][pijiu][lanqiu][juhua][hecai][haobang][caidao][baojin][chi][dan][kulou][shuai][shouqiang][yangtuo][youling]
代码

请输入代码：

确认
图片

    暂无评论内容