低通滤波器和高通滤波器

应用于图像低通滤波器和高通滤波器的实现

需要用到傅里叶变换


#include <opencv2/opencv.hpp>
#include <Eigen>
#include <iostream>
#include <vector>
#include <cmath>
#include <complex>#define M_PI       3.14159265358979323846   // pi
// 对数幅度缩放
Eigen::MatrixXd logAmplitudeSpectrum(const Eigen::MatrixXd& spectrum) {return (spectrum.array() + 1).log();
}//  乘幂尺度变换
Eigen::MatrixXd powerLawScaling(const Eigen::MatrixXd& spectrum, double gamma) {return spectrum.array().pow(gamma);
}// 归一化 to [0, 1]
Eigen::MatrixXd normalize(const Eigen::MatrixXd& spectrum) {double minVal = spectrum.minCoeff();double maxVal = spectrum.maxCoeff();return (spectrum.array() - minVal) / (maxVal - minVal);
}// 增强频谱显示
Eigen::MatrixXd enhanceSpectrumDisplay(const Eigen::MatrixXd& spectrum, double gamma = 1) {Eigen::MatrixXd logSpectrum = logAmplitudeSpectrum(spectrum);Eigen::MatrixXd powerScaledSpectrum = powerLawScaling(logSpectrum, gamma);return normalize(powerScaledSpectrum);
}// 从复数矩阵中获取幅度谱和相位谱
void getAmplitudeAndPhaseSpectra(const Eigen::MatrixXcd& data, Eigen::MatrixXd& amplitude, Eigen::MatrixXd& phase) {amplitude = data.array().abs().matrix();phase = data.array().arg().matrix();
}// 从幅度谱和相位谱重构复数矩阵
void reconstructFromAmplitudeAndPhase(const Eigen::MatrixXd& amplitude,const Eigen::MatrixXd& phase,Eigen::MatrixXcd& data)
{data = (amplitude.array() * (phase.array().cos() + std::complex<double>(0, 1) * phase.array().sin())).matrix();}// 1:**振幅谱低频移动到中心(频率平移)**:方便操作,利用象限对称互换 fft之后
Eigen::MatrixXd fftShift(const Eigen::MatrixXd &F)
{int M = F.rows();int N = F.cols();Eigen::MatrixXd F_shifted(M, N);int mid_M = M >> 1;int mid_N = N >> 1;// 交换第一象限和第三象限F_shifted.block(0, 0, mid_M, mid_N) = F.block(mid_M, mid_N, mid_M, mid_N);F_shifted.block(mid_M, mid_N, mid_M, mid_N) = F.block(0, 0, mid_M, mid_N);// 交换第二象限和第四象限F_shifted.block(0, mid_N, mid_M, mid_N) = F.block(mid_M, 0, mid_M, mid_N);F_shifted.block(mid_M, 0, mid_M, mid_N) = F.block(0, mid_N, mid_M, mid_N);return F_shifted;
}// 2:振幅谱低频移动到中心  **图像进行-1幂操作**:然后经过fft变换后,低频会在振幅谱中间 fft之前
Eigen::MatrixXd imageShift(const Eigen::MatrixXd &image)
{int M = image.rows();int N = image.cols();Eigen::MatrixXd F_shifted(M, N);// 通过乘以 (-1)^(u+v) 来平移频率for (int u = 0; u < M; ++u){for (int v = 0; v < N; ++v){//(u + v) & 1 通过位的判断末尾如果是1 为奇数,为0,为偶数。 -1的奇次幂还是-1,偶次幂为1.F_shifted(u, v) = image(u, v) * ((u + v) & 1 ? -1 : 1);}}return F_shifted;
}Eigen::MatrixXd tMatrixXd(const cv::Mat &img)
{Eigen::MatrixXd image(img.rows, img.cols);for (int i = 0; i < img.rows; ++i){for (int j = 0; j < img.cols; ++j){image(i, j) = img.at<float>(i, j);}}return image;
}
cv::Mat tMat(const Eigen::MatrixXd &img)
{cv::Mat image(img.rows(), img.cols(), CV_32FC1);for (int i = 0; i < img.rows(); ++i){for (int j = 0; j < img.cols(); ++j){image.at<float>(i, j) = img(i, j);}}return image;
}//inv =true 为逆变换,false为 正变换
//FFT 使用Eigen库中的向量表示,方便二维计算
Eigen::VectorXcd FFT(const Eigen::VectorXcd& y, bool inv = false)
{Eigen::VectorXcd x = y;//数据的大小int N = x.size();// 按位反转for (int i = 1, j = 0; i < N; i++){//bit = N/2int bit = N >> 1;/** 我们正在查看j的二进制表示中的每一位。从左到右,我们检查每一位是否为1。对于每一个为1的位,我们将其反转为0。当我们遇到第一个为0的位时,循环终止。*/for (; j & bit; bit >>= 1){j ^= bit;}//进行异或运算(XOR运算的工作原理是:当两个比较的位相同时,结果是0;当两个比较的位不同时,结果是1。)j ^= bit;//当 i >j 表示前面已经替换了位置, i==j,表示位置不用变if (i < j){std::swap(x[i], x[j]);}}// 预先计算旋转因子Eigen::VectorXcd w(N >> 1);//逆变换 = 1;傅里叶变换  = -1;double imag_i = inv ? 1.0 : -1.0;std::complex<double> tempW = std::exp(std::complex<double>(0, imag_i * 2 * M_PI / N));w[0] = 1;for (int i = 1; i < (N >> 1); ++i){// w[i] = std::polar(1.0, imag_i * 2 * M_PI * i / N);//w[i] = std::pow(tempW, i);w[i] = w[i - 1] * tempW;}// 迭代FFTfor (int len = 2; len <= N; len <<= 1){int halfLen = len >> 1;int step = N / len;for (int i = 0; i < N; i += len){for (int j = 0; j < halfLen; ++j){std::complex<double> u = x[i + j];std::complex<double> v = x[i + j + halfLen] * w[j * step];x[i + j] = u + v;x[i + j + halfLen] = u - v;}}}//使用逆变换时if (inv){for (std::complex<double>& a : x) {a /= N;}}return x;
}//检查二维数据大小是否为2的幂数,不是则填充为0到2的幂数大小
Eigen::MatrixXd padToPowerOfTwo(const Eigen::MatrixXd& matrix) {int rows = matrix.rows();int cols = matrix.cols();// 确保输入的大小是2的幂int newRows = 1, newCols = 1;while (newRows < rows){newRows <<= 1;}while (newCols < cols){newCols <<= 1;}Eigen::MatrixXd paddedMatrix = Eigen::MatrixXd::Zero(newRows, newCols);paddedMatrix.block(0, 0, rows, cols) = matrix;return paddedMatrix;
}
// 离散傅里叶变换 -  二维
Eigen::MatrixXcd FFT2D(const Eigen::MatrixXcd& image, bool inv = false) {int rows = image.rows();int cols = image.cols();Eigen::MatrixXcd result(rows, cols);for (int i = 0; i < rows; i++) {Eigen::VectorXcd temp = image.row(i).transpose();result.row(i) = FFT(temp, inv);}for (int j = 0; j < cols; j++) {Eigen::VectorXcd temp = result.col(j);result.col(j) = FFT(temp, inv);}return result;
}//创建高频  radius越小,越减少高频的部分,越大,越还原图像
Eigen::MatrixXd highFrequency (const Eigen::MatrixXd & data,int radius)
{int rows = data.rows();int cols = data.cols();// 创建高通滤波器(圆形掩码)Eigen::MatrixXd mask = Eigen::MatrixXd::Ones(rows, cols);//中心坐标int centerRow = rows >> 1;int centerCol = cols >> 1;//找到以 radius 大小的矩形范围内//左边位置int left = centerCol - radius;//超过边界 为0left = left > 0 ? left : 0;//右边位置int right = centerCol + radius;//超过边界 为0right = right < cols ? right : cols;//上边位置int top = centerRow - radius;//超过边界 为0top = top > 0 ? top : 0;//下边位置int down = centerRow + radius;//超过边界 为0down = down < rows ? down : rows;//在正矩形内画最大的圆for (int i = top; i < down; ++i){for (int j = left; j < right; ++j){//在图像中心画圆,半径不能超过double distance = std::sqrt(std::pow(i - centerRow, 2) + std::pow(j - centerCol, 2));if (distance <= radius){mask(i, j) = 0.0;}}}return mask;}
//创建高通滤波器 -
//简单的说,就是靠近频谱图中心的低频部分给舍弃掉,远离频谱图中心的高频部分保留。通常会保留物体的边界。
Eigen::MatrixXd  highPassFilter(const Eigen::MatrixXd & image, int radius = 0)
{int rows = image.rows();int cols = image.cols();//大小变为2的幂数Eigen::MatrixXd data = padToPowerOfTwo(image);// 计算傅里叶变换Eigen::MatrixXcd transformed = FFT2D(data);//获取幅度谱和相位谱Eigen::MatrixXd amplitude, phase;getAmplitudeAndPhaseSpectra(transformed, amplitude, phase);//振幅谱移动到中心(频率平移)amplitude = fftShift(amplitude);/////增强振幅,用于观测 --  实际运算注释掉Eigen::MatrixXd amplitude1 = enhanceSpectrumDisplay(amplitude,1);cv::Mat highP = tMat(amplitude1);cv::imshow("highPassFilter",highP);/////根据radius创建高频掩码Eigen::MatrixXd mask = highFrequency(amplitude, radius);amplitude = amplitude.array() * mask.array();//振幅谱移动到中心(频率平移)反转换amplitude = fftShift(amplitude);// 从幅度谱和相位谱重构复数矩阵reconstructFromAmplitudeAndPhase(amplitude, phase, transformed);// 计算逆变换Eigen::MatrixXcd reconstructed = FFT2D(transformed, true);return reconstructed.real().block(0, 0, rows, cols);
}
//创建低频  radius越小,越还原图像,越大,减少低频的部分,
Eigen::MatrixXd lowFrequency(const Eigen::MatrixXd & data,int radius)
{int rows = data.rows();int cols = data.cols();// 创建低通滤波器(圆形掩码)Eigen::MatrixXd mask = Eigen::MatrixXd::Zero(rows, cols);//中心坐标int centerRow = rows >> 1;int centerCol = cols >> 1;//找到以 radius 大小的矩形范围内//左边位置int left = centerCol - radius;//超过边界 为0left = left > 0 ? left : 0;//右边位置int right = centerCol + radius;//超过边界 为0right = right < cols ? right : cols;//上边位置int top = centerRow - radius;//超过边界 为0top = top > 0 ? top : 0;//下边位置int down = centerRow + radius;//超过边界 为0down = down < rows ? down : rows;//在正矩形内画最大的圆for (int i = top; i < down; ++i){for (int j = left; j < right; ++j){//在图像中心画圆,半径不能超过double distance = std::sqrt(std::pow(i - centerRow, 2) + std::pow(j - centerCol, 2));if (distance <= radius){mask(i, j) = 1.0;}}}return mask;}
//创建低通滤波器 -
//简单的说,就是靠近频谱图中心的低频部分给保留,远离频谱图中心的高频部分给去除掉。但是这会影响图像的清晰度。
Eigen::MatrixXd  lowPassFilter(const Eigen::MatrixXd & image, int radius = 0)
{int rows = image.rows();int cols = image.cols();//大小变为2的幂数Eigen::MatrixXd data = padToPowerOfTwo(image);// 计算傅里叶变换Eigen::MatrixXcd transformed = FFT2D(data);//获取幅度谱和相位谱Eigen::MatrixXd amplitude, phase;getAmplitudeAndPhaseSpectra(transformed, amplitude, phase);//振幅谱移动到中心(频率平移)amplitude = fftShift(amplitude);///
//    //增强振幅,用于观测 --  实际运算注释掉
//    Eigen::MatrixXd amplitude1 = enhanceSpectrumDisplay(amplitude,1);
//    cv::Mat highP = tMat(amplitude1);
//    cv::imshow("lowPassFilter",highP);/////根据radius创建高频掩码Eigen::MatrixXd mask = lowFrequency(amplitude, radius);amplitude = amplitude.array() * mask.array();//振幅谱移动到中心(频率平移)反转换amplitude = fftShift(amplitude);// 从幅度谱和相位谱重构复数矩阵reconstructFromAmplitudeAndPhase(amplitude, phase, transformed);// 计算逆变换Eigen::MatrixXcd reconstructed = FFT2D(transformed, true);return reconstructed.real().block(0, 0, rows, cols);
}int main()
{cv::Mat img = cv::imread("193560523230866.png");if (img.empty()){std::cout << "请确定是否输入正确的图像文件" << std::endl;}cv::Mat gray;cvtColor(img, gray, cv::COLOR_BGR2GRAY);//图像转换CV_32F储存gray.convertTo(gray, CV_32F, 1 / 255.0, 0);//图像太大,用直接计算 耗时太长,缩小比例//resize(gray, gray, cv::Size(80, 80));//Mat  转 MatrixXdEigen::MatrixXd image = tMatrixXd(gray);// 记录开始时间auto start = std::chrono::high_resolution_clock::now();//低频Eigen::MatrixXd low = lowPassFilter(image,50);//高频Eigen::MatrixXd high = highPassFilter(image, 50);// 记录结束时间auto stop = std::chrono::high_resolution_clock::now();// 计算持续时间auto duration = std::chrono::duration_cast<std::chrono::microseconds>(stop - start);qDebug() << "代码运行时长: " << duration.count() << " 微秒" ;cv::Mat lowP = tMat(low);//+0.5 增加显示效果high = high.array() + 0.5;cv::Mat highP = tMat(high);cv::imshow("lowP",lowP);cv::imshow("highP",highP);return 0;
}

在这里插入图片描述

在这里插入图片描述

  • 振幅增强
    在这里插入图片描述
    在这里插入图片描述

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/54083.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

QT5.12.12通过ODBC连接到GBase 8s数据库(CentOS)

本示例使用的环境如下&#xff1a; 硬件平台&#xff1a;x86_64&#xff08;amd64&#xff09;操作系统&#xff1a;CentOS 7.8 2003数据库版本&#xff08;含CSDK&#xff09;&#xff1a;GBase 8s V8.8 3.0.0_1 为什么使用QT 5.12.10&#xff1f;该版本包含QODBC。 1&#…

ES6中promise的使用

ES6中promise的使用 本文目录 ES6中promise的使用基础介绍箭头函数function函数状态 原型方法Promise.prototype.then()Promise.prototype.catch() 静态方法Promise.all()Promise.race()Promise.any() 链式回调 基础介绍 官网&#xff1a;https://promisesaplus.com/ window.…

最新docker多系统安装技术

在Ubuntu操作系统中安装Docker 在Ubuntu操作系统中安装Docker的步骤如下。 1&#xff0e;卸载旧版本Docker 卸载旧版本Docker的命令如下&#xff1a; $ sudo apt-get remove docker docker-engine docker.io 2&#xff0e;使用脚本自动安装 在测试或开发环境中&#xff0…

STM32 进不了main 函数

1. 我用的是STM32L151C8T6 的芯片&#xff0c;在github 上找了个别人的例程&#xff0c;拿来当模板改&#xff0c;由于他用的是HSE 外部晶振&#xff0c;我用的是内部晶振HSI&#xff0c;所以需要改系统时钟&#xff0c;改完后debug&#xff0c; 一直进不了main 函数&#xff0…

PHP“牵手”拼多多商品详情数据获取方法,拼多多API接口批量获取拼多多商品详情数据说明

拼多多商品详情接口 API 是开放平台提供的一种 API 接口&#xff0c;它可以帮助开发者获取拼多多商品的详细信息&#xff0c;包括商品的标题、描述、图片等信息。在拼多多电商平台的开发中&#xff0c;拼多多详情接口 API 是非常常用的 API&#xff0c;因此本文将详细介绍拼多多…

【C++】C++ 引用详解 ⑤ ( 函数 “ 引用类型返回值 “ 当左值被赋值 )

文章目录 一、函数返回值不能是 " 局部变量 " 的引用或指针1、函数返回值常用用法2、分析函数 " 普通返回值 " 做左值的情况3、分析函数 " 引用返回值 " 做左值的情况 函数返回值 能作为 左值 , 是很重要的概念 , 这是实现 " 链式编程 &quo…

淘宝API技术解析,实现关键词搜索淘宝商品(商品详情接口等)

淘宝提供了开放平台接口&#xff08;API&#xff09;来实现按图搜索淘宝商品的功能。您可以通过以下步骤来实现&#xff1a; 获取开放平台的访问权限&#xff1a;首先&#xff0c;您需要在淘宝开放平台创建一个应用&#xff0c;获取访问淘宝API的权限。具体的申请步骤和要求可以…

LabVIEW开发灭火器机器人

LabVIEW开发灭火器机器人 如今&#xff0c;自主机器人在行业中有着巨大的需求。这是因为它们根据不同情况的适应性。由于消防员很难进入高风险区域&#xff0c;自主机器人出现了。该机器人具有自行检测火灾的能力&#xff0c;并通过自己的决定穿越路径。 由于消防安全是主要问…

java八股文面试[java基础]——如何实现不可变的类

知识来源&#xff1a; 【23版面试突击】如何实现不可变的类&#xff1f;_哔哩哔哩_bilibili 【2023年面试】怎样声明一个类不会被继承&#xff0c;什么场景下会用_哔哩哔哩_bilibili

JVM 之字节码(.class)文件

本文中的内容参考B站尚硅谷宋红康JVM全套教程 你将获得&#xff1a; 1、掌握字节码文件的结构 2、掌握Java源代码如何在JVM中执行 3、掌握一些虚拟机指令 4、回答一些面试题 课程介绍 通过几个面试题初始字节码文件为什么学习class字节码文件什么是class字节码文件分析c…

2022年03月 C/C++(四级)真题解析#中国电子学会#全国青少年软件编程等级考试

第1题&#xff1a;拦截导弹 某国为了防御敌国的导弹袭击&#xff0c; 发展出一种导弹拦截系统。 但是这种导弹拦截系统有一个缺陷&#xff1a; 虽然它的第一发炮弹能够到达任意的高度&#xff0c;但是以后每一发炮弹都不能高于前一发的高度。 某天&#xff0c; 雷达捕捉到敌国的…

Vue3.0极速入门- 目录和文件说明

目录结构 以下文件均为npm create helloworld自动生成的文件目录结构 目录截图 目录说明 目录/文件说明node_modulesnpm 加载的项目依赖模块src这里是我们要开发的目录&#xff0c;基本上要做的事情都在这个目录里assets放置一些图片&#xff0c;如logo等。componentsvue组件…

SFM structure from motion

struction就是空间三维点的位置 motion 就是相机每帧的位移 https://www.youtube.com/watch?vUhkb8Zq-dnM&listPL2zRqk16wsdoYzrWStffqBAoUY8XdvatV&index9

西部AI小镇-构建自主虚拟世界

背景 未来曜文有接入市场上所有面向chatGPT开发的应用&#xff0c;例如开源聊天组件&#xff0c;西部小镇等 内容介绍 生成代理起床&#xff0c;做早餐&#xff0c;然后去上班&#xff1b;艺术家作画&#xff0c;作家写作&#xff1b;他们形成意见、互相关注并发起对话&…

Linux线程 --- 生产者消费者模型(C语言)

在学习完线程相关的概念之后&#xff0c;本节来认识一下Linux多线程相关的一个重要模型----“ 生产者消费者模型” 本文参考&#xff1a; Linux多线程生产者与消费者_红娃子的博客-CSDN博客 Linux多线程——生产者消费者模型_linux多线程生产者与消费者_两片空白的博客-CSDN博客…

基于Java+SpringBoot+Vue前后端分离党员教育和管理系统设计和实现

博主介绍&#xff1a;✌全网粉丝30W,csdn特邀作者、博客专家、CSDN新星计划导师、Java领域优质创作者,博客之星、掘金/华为云/阿里云/InfoQ等平台优质作者、专注于Java技术领域和毕业项目实战✌ &#x1f345;文末获取源码联系&#x1f345; &#x1f447;&#x1f3fb; 精彩专…

ServiceManager接收APP的跨进程Binder通信流程分析

现在一起来分析Server端接收&#xff08;来自APP端&#xff09;Binder数据的整个过程&#xff0c;还是以ServiceManager这个Server为例进行分析,这是一个至下而上的分析过程。 在分析之前先思考ServiceManager是什么&#xff1f;它其实是一个独立的进程&#xff0c;由init解析i…

银河麒麟服务器、centos7服务器一键卸载mysql脚本

脚本 # 查看mysql相关的rpm包写到rmsql.sh文件中 rpm -aq | grep -i mysql >rmsql.sh # 修改文件为卸载mysql的脚本文件 sed -i -e s/^/yum remove -y / rmsql.sh # 修改文本权限 chmod 777 rmsql.sh # 全盘查找mysql相关文件&#xff0c;写到my.sh脚本中 find / -name mysq…

git及GitHub的使用

文章目录 git在本地仓库的使用github使用创建仓库https协议连接(不推荐&#xff0c;现在用起来比较麻烦)ssh连接&#xff08;推荐&#xff09;git分支操作冲突处理忽略文件 git在本地仓库的使用 1.在目标目录下右键打开git bash here 2.创建用户名和邮箱(注&#xff1a; 下载完…

框架(Git基础详解及Git在idea中集成步骤)

目录 基础&#xff1a; idea集成Git并添加项目到git仓库 1.idea集成git&#xff0c;集成.git.exe文件 2.初始化本地Git仓库项目 3. 将工作区代码添加到暂存区 4.将暂存区代码添加到本地仓库 5.Git本地库操作 Idea集成Gitee并提交代码到第三方库 1.setting里搜索gitee 2.添…