C++版本的OpenCV实现二维图像的卷积定理(通过傅里叶变换实现二维图像的卷积过程,附代码!!)

2023-09-15 10:21:22

前言

工作中用到许多卷积过程,需要转成C++代码的实现,使用OpenCV库自带的二维卷积过程所耗费的时间比较久,为了提升代码的运行效率可以考虑使用卷积定理实现二维图像的卷积过程。

一、卷积定理简单介绍

卷积定理是傅立叶变换满足的一个重要性质。卷积定理指出,函数卷积的傅立叶变换是函数傅立叶变换的乘积。具体分为时域卷积定理和频域卷积定理,时域卷积定理即时域内的卷积对应频域内的乘积;频域卷积定理即频域内的卷积对应时域内的乘积,两者具有对偶关系。

二、不同卷积过程对应的傅里叶变换过程

1、“Same”卷积

假设原始图像的大小为mxm,卷积核的大小为nxn。此时进行same卷积,则卷积后生成的图像大小为mxm,此时进行卷积是需要对原始图像进行padding操作,需要对原始图像周围进行n-1个补零操作,此时被卷积图像大小padding为(m+n-1)x(m+n-1),然后进行卷积操作,进行傅里叶变化是需要先将被卷积图像以及卷积核的大小都padding为(m+n-1)x(m+n-1),因为要进行傅里叶变换后的乘积操作,因此被卷积图像以及卷积核的大小需要相等。经过卷积定理的操作之后,将生成的(m+n-1)x(m+n-1)大小的图像按照padding的逆操作进行裁剪,得到mxm的图像即为same卷积得到的卷积图像。
Python代码验证

import numpy as np
from scipy import signal

# 原始图像 f(x)
gray = np.uint16(np.random.randint(100, size=(7, 7)))
# 卷积核 g(x)
kenel = np.ones((3, 3))/9
# ----- Conv = f(x)*g(x) ----- #
# f(x)*g(x)
Conv = signal.convolve2d(gray, kenel, mode='same') # 使用full卷积类型,得到(M+N-1)X(M+N-1)大小

#--------- ifft{ F(f(x))·F(g(x)) } ---------#
# 傅里叶变换前进行 padding 填充,图像和卷积核都补零到 (M+N-1)x(M+N-1)大小。
img_pad = np.pad(gray, ((1, 1), (1, 1)), 'constant')
kenel_pad = np.pad(kenel, ((3, 3),(3, 3)), 'constant')

# F(f(x))
img_fft = np.fft.fftshift(np.fft.fft2(img_pad))
# F(g(x))
kenel_fft = np.fft.fftshift(np.fft.fft2(kenel_pad))
# ifft( F(f(x))·F(g(x)) )
FFT = np.fft.ifftshift(np.fft.ifft2(np.fft.fftshift(img_fft*kenel_fft)))

#--------- 打印结果 ---------#
print(" f(x) ↓")
print(gray)
print(" g(x) ↓")
print(kenel)

print("\n\n f(x)*g(x) ↓")
print(np.uint8(Conv))
print("\n ifft[F·G] ↓")
print(np.uint8(np.abs(FFT)))

结果展示
上述代码的结果如下图所示,可见卷积定理得到的结果经过裁剪后的红框中的内容与same卷积得到的结果一致。
在这里插入图片描述

2、“Full”卷积

假设原始图像的大小为mxm,卷积核的大小为nxn。此时进行full卷积,则卷积后生成的图像大小为(m+n-1)x(m+n-1),此时进行卷积是需要对原始图像进行padding操作,需要对原始图像周围进行2n-2个补零操作,此时被卷积图像大小padding为(m+2n-2)x(m+2n-12),然后进行卷积操作,进行傅里叶变化是需要先将被卷积图像以及卷积核的大小都padding为(m+n-1)x(m+n-1),因为要进行傅里叶变换后的乘积操作,因此被卷积图像以及卷积核的大小需要相等。经过卷积定理的操作之后,得到(m+n-1)x(m+n-1)的图像即为full卷积得到的卷积图像。
Python代码验证

import numpy as np
from scipy import signal

# 原始图像 f(x)
gray = np.uint16(np.random.randint(100, size=(9, 9)))
# 卷积核 g(x)
kenel = np.ones((5, 5))/9
# ----- Conv = f(x)*g(x) ----- #
# f(x)*g(x)
Conv = signal.convolve2d(gray, kenel, mode='full') # 使用full卷积类型,得到(M+N-1)X(M+N-1)大小


#--------- ifft{ F(f(x))·F(g(x)) } ---------#
# 傅里叶变换前进行 padding 填充,图像和卷积核都补零到 (M+N-1)x(M+N-1)大小。
img_pad = np.pad(gray, ((2, 2), (2, 2)), 'constant')
kenel_pad = np.pad(kenel, ((4, 4), (4, 4)), 'constant')

# F(f(x))
img_fft = np.fft.fftshift(np.fft.fft2(img_pad))
# F(g(x))
kenel_fft = np.fft.fftshift(np.fft.fft2(kenel_pad))
# ifft( F(f(x))·F(g(x)) )
FFT = np.fft.ifftshift(np.fft.ifft2(np.fft.fftshift(img_fft*kenel_fft)))

#--------- 打印结果 ---------#
print(" f(x) ↓")
print(gray)
print(" g(x) ↓")
print(kenel)

print("\n\n f(x)*g(x) ↓")
print(np.uint8(Conv))
print("\n ifft[F·G] ↓")
print(np.uint8(np.abs(FFT)))

结果展示
上述代码的结果如下图所示,可见卷积定理得到的结果与full卷积得到的结果一致。
在这里插入图片描述

3、“Valid”卷积

假设原始图像的大小为mxm,卷积核的大小为nxn。此时进行valid卷积,则卷积后生成的图像大小为(m-n+1)x(m-n+1),此时不需要对被卷积图像进行padding操作,直接进行卷积操作,进行傅里叶变化是需要先将卷积核的大小padding为mxm,因为要进行傅里叶变换后的乘积操作,因此被卷积图像以及卷积核的大小需要相等。经过卷积定理的操作之后,将生成的mxm的图像进行裁剪,得到(m-n+1)x(m-n+1)的图像即为valid卷积得到的卷积图像。
Python代码验证

import numpy as np
from scipy import signal

# 原始图像 f(x)
gray = np.uint16(np.random.randint(100, size=(9, 9)))
# 卷积核 g(x)
kenel = np.ones((5, 5))/9
# ----- Conv = f(x)*g(x) ----- #
# f(x)*g(x)
Conv = signal.convolve2d(gray, kenel, mode='valid') # 使用full卷积类型,得到(M+N-1)X(M+N-1)大小

#--------- ifft{ F(f(x))·F(g(x)) } ---------#
# 傅里叶变换前进行 padding 填充,图像和卷积核都补零到 (M+N-1)x(M+N-1)大小。
img_pad = np.pad(gray, ((0, 0), (0, 0)), 'constant')
kenel_pad = np.pad(kenel, ((2, 2), (2, 2)), 'constant')

# F(f(x))
img_fft = np.fft.fftshift(np.fft.fft2(img_pad))
# F(g(x))
kenel_fft = np.fft.fftshift(np.fft.fft2(kenel_pad))
# ifft( F(f(x))·F(g(x)) )
FFT = np.fft.ifftshift(np.fft.ifft2(np.fft.fftshift(img_fft*kenel_fft)))

#--------- 打印结果 ---------#
print(" f(x) ↓")
print(gray)
print(" g(x) ↓")
print(kenel)

print("\n\n f(x)*g(x) ↓")
print(np.uint8(Conv))
print("\n ifft[F·G] ↓")
print(np.uint8(np.abs(FFT)))

结果展示
上述代码的结果如下图所示,可见卷积定理得到的结果经过裁剪后的红框中的内容与valid卷积得到的结果一致。
在这里插入图片描述

三、基于OpenCV库实现的二维图像卷积定理

//本代码实现的是使用193x193大小的卷积核对193x193大小的被卷积图像进行卷积操作
Conv2Kernel(kernelImage, inputImage, outputImage, n, n);
cv::Mat matPaded1;
cv::Mat matPaded2;
cv::Mat kernelPaded;
cv::Mat fftKernel;
cv::Mat fftMat1;
cv::Mat fftMat2;
cv::Mat result1;
cv::Mat result2;

cv::copyMakeBorder(matKernel, kernelPaded, 192, 215, 192, 215, cv::BORDER_CONSTANT, cv::Scalar(0.f));
cv::copyMakeBorder(matConv1, matPaded1, 192, 215, 192, 215, cv::BORDER_CONSTANT, cv::Scalar(0.f));

cv::Mat planes1[] = { cv::Mat_<float>(kernelPaded),cv::Mat::zeros(kernelPaded.size(),CV_32F)
};
cv::Mat planes2[] = { cv::Mat_<float>(matPaded1),cv::Mat::zeros(matPaded1.size(),CV_32F)
};
cv::merge(planes1, 2, fftKernel);
cv::merge(planes2, 2, fftMat1);

cv::dft(fftKernel, fftKernel, cv::DFT_COMPLEX_OUTPUT);
cv::dft(fftMat1, fftMat1, cv::DFT_COMPLEX_OUTPUT);

fftshift(fftKernel);
fftshift(fftMat1);
cv::Mat fftMultiplication1;
cv::mulSpectrums(fftKernel, fftMat1, fftMultiplication1, cv::DFT_ROWS);
cv::idft(fftMultiplication1, result1, cv::DFT_INVERSE + cv::DFT_SCALE + cv::DFT_COMPLEX_OUTPUT);

ifftshift(result1);
cv::split(result1, planes1);
cv::magnitude(planes1[0], planes1[1], planes1[0]);

cv::Mat matConv3 = planes1[0](cv::Rect(192, 192, 193, 193)).clone();

代码详解
1、对被卷积图像和卷积核执行padding操作
由于使用的是same卷积,按照上文same卷积的过程对其进行padding操作,图像上下左右四个方向都补上192个0。下述代码右边和下面之所以补零数不是192,是因为在进行傅里叶变换时,特定长度的矩阵速度更快,因此补了更多的0。但是我在实验过程中并没有发现速度有变快
2、对padding后的图像矩阵进行傅里叶变换
对padding后的卷积核以及被卷积图像执行傅里叶变换,然后进行fftshift操作。
3、对傅里叶变换后的图像矩阵进行点乘操作
4、对点乘后的举证进行傅里叶逆变换
对点乘后的矩阵执行傅里叶逆变换,然后进行ifftshift操作。
5、截取合适位置的图像,具体位置根据padding过程确定,得到的图像矩阵即为卷积后的图像矩阵。

四、基于FFTW库实现的二维图像卷积定理

                fftw_complex* pImgOut = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * 600 * 600);
                fftw_plan pRef;
                pRef = fftw_plan_dft_2d(600, 600, pImgIn, pImgOut, FFTW_FORWARD, FFTW_ESTIMATE);

                fftw_execute(pRef);
                fftw_destroy_plan(pRef);
                fftw_free(pImgIn);
                fftw_free(pImgOut);

五、总结与讨论

经过不懈努力终于使用C++实现了卷积定理,想要进一步提升卷积过程的速度,优化代码性能。但是天不遂人愿,在我的任务中使用卷积定义实现卷积过程和使用OpenCV库中的卷积操作接口速度相差无几,几乎没有看到在计算速度上有提升,忙了半天白忙活了,难受~。分析原因可能是由于本次任务的卷积过程比较特殊,卷积核、被卷积图像的大小一致,导致进行padding操作时,padding后的图像几乎是原图像的9倍大小,造成傅里叶变换以及傅里叶逆变换的过程速度变慢,影响了整体流程的速度。

更多推荐

python基础

一,什么是pythonPython是一种高级、通用且解释型的编程语言,由GuidovanRossum于1991年首次发布。它具有简洁的语法、清晰的代码结构和强大的功能,被广泛应用于各种领域,包括软件开发、数据分析、人工智能、网络编程等。以下是Python的一些特点和优势:简单易学:Python具有直观、简洁的语法,易于

GDPU 数据结构 天码行空2

实验内容用顺序表实现病历信息的管理与查询功能。具体要求如下:利用教材中定义顺序表类型存储病人病历信息(病历号,姓名,症状);要求使用头文件。设计顺序表定位查找算法,写成一个函数,完成的功能为:在线性表L中查找数据元素x,如果存在则返回线性表中和x值相等的第1个数据元素的序号;如果不存在,则返回-1。函数定义为intLi

数据结构之-----二叉树

目录本章内容如下:1:树的相关概念与结构2:二叉树的概念与结构3:二叉树的链式结构与实现文章正式开始,让我们一起学习树吧!!一:树的概念树是一种非线性结构,与我们前面所学的顺序表与链表不同,数据元素的对应是1对多的关系,只有一个根结点,且除了根节点其它的结点有且仅有1个前驱结点(父结点)。我们可以将一棵树看作由很多个结

计算机和编程语言初见

学习程序设计的目的是什么呢?不一定要做出一个软件或系统出来,更重要的是理解计算机是如何工作的以及它的长处和短处。计算机本身是无意识的,因此我们要求它为我们做事时:应该将步骤细化、“直”化(规律化);其实计算机什么也不会,我们必须手把手地教他一步一步的做。而计算机的某个优点也正是如此——听话,你叫它往东它绝不往西。然后我

Stellar Toolkit for MySQL 9.0 Crack 3in1

面向数据库管理员的MySQL工具包StellarToolkitforMySQL是一款三合一软件套件,用于修复损坏的MySQL和MariaDB数据库、从MySQL数据库的InnoDB和MyISAM表恢复数据以及分析MySQL数据库日志文件。该软件还可以以最高的安全性和完整性相互转换MySQL/MariaDB、MSSQL(

【跟晓月学数据库】基于book库的mysql进阶实战

前言上篇文章中,我们已经导入了book库,如果你还没有导入book库,参考:【跟晓月学数据库】使用MySQLdump对数据导入导出这篇文章,主要是基于book库的操作,希望对你有用。🏠个人主页:我是沐风晓月🧑个人简介:大家好,我是沐风晓月,阿里云社区专家博主😉😉💕座右铭:先努力成长自己,再帮助更多的人,一起加

运营商大数据精准营销获客?

多年来,大数据运营商一直致力于为企业提供互联网大数据精准营销的新项目,并以确保自身信息安全为前提。例如,如果移动用户查看了任何网站,在网页上搜索了任何关键词,登录了应用程序,给任何人打了电话,以及隶属地区、性別,所有这些都由运营商存储,那么企业可以提供需求,运营商可以根据客户行为找到准确的意向客户。高质量的新客户可以基

定义爬虫规则和数据存储

定义爬虫规则是指确定爬虫程序应该如何访问和提取网页数据的规则。这些规则包括确定要爬取的网页的URL、确定要提取的数据类型和位置、确定爬取的深度和频率等。爬虫规则通常由以下几个方面组成:起始URL:确定爬虫程序开始爬取的网页URL。URL过滤规则:确定哪些URL应该被爬取,哪些URL应该被忽略。可以使用正则表达式或其他方

支持向量机(SVM)案例分析

支持向量机(supportvectormachines,SVM)是一种二分类模型,所谓二分类模型是指比如有很多特征(自变量X)对另外一个标签项(因变量Y)的分类作用关系,比如当前有很多特征,包括身高、年龄、学历、收入、教育年限等共5项,因变量为‘是否吸烟’,‘是否吸烟’仅包括两项,吸烟和不吸烟。那么该5个特征项对于‘是

【C++基于多设计模式下的同步&异步日志系统】

文章目录@[toc]1:peach:项目介绍:peach:2:peach:开发环境:peach:3:peach:核心技术:peach:4:peach:环境搭建:peach:5:peach:日志系统介绍:peach:5.1:apple:为什么需要日志系统?:apple:5.2:apple:日志系统技术实现:apple:5

Java 消息策略的实现 - Kafak 是怎么设计的

这个也是开放讨论题,主要讨论下Kafka在消息中是如何进行实现的。1_cCyPNzf95ygMFUgsrleHtw976×50621.4KB总结这个题目的开发性太强了。Kafka可以用的地方非常多,我经历过的项目有Kafka用在消息处理策略上的。这个主要是IoT项目,因为这个项目需要对温度传感器采集获得数据。当我们有多

热文推荐