本文从一名工程师的角度讲解了数字图像处理领域的重要概念——卷积。本文力求使用最少的数学公式将卷积的物理含义以及计算方法讲解清楚。其中第四节和第五节图解卷积过程,第六节给出了卷积的代码实现。相信工程师朋友们可以通过这一篇文章轻松理解卷积。
(一)线性时不变系统
如果一个系统具备齐次性,可加性和时不变性,则该系统为线性时不变系统。
齐次性:如果一个系统对输入信号
可加性:如果一个系统对输入信号
时不变性:如果一个系统对输入信号
(二)单位脉冲信号
如果一个信号
信号
任何离散信号都可以表示为单位脉冲信号的线性组合。如下图所示,信号x[n]可以表示为一组脉冲信号的线性组合:
(三)单位脉冲响应
如果已知一个线性时不变系统的单位脉冲响应(输入信号为单位脉冲信号
(1)因为任何离散信号均可以表示成单位脉冲信号的线性组合。
(2)由线性时不变系统的齐次性和可加性可知,线性时不变系统的输出为:
(3)由线性时不变系统的时不变性可知:
(4)已知线性时不变系统的单位脉冲响应为
上式即为卷积,即已知线性时不变系统的单位脉冲响应
(四)卷积计算(视角一)
该视角从输入信号出发,对于输入信号中的每一个元素
当输入信号的所有元素均完成该操作后,将所有的结果信号相加,即可得到
视角一的每次迭代,输入信号中的一个元素就完成了它的使命,输入信号
(五)卷积计算(视角二)
该视角从输出信号y[n]出发,根据视角一中的分析,可以观察到输出信号的每个元素的计算方法如下图所示:
该视角每次迭代完成输出信号中的一个元素的计算,要注意该方法中信号
(六)卷积计算代码实现
本节从文中的视角一和视角二两个角度实现卷积的计算。卷积的计算类的代码粘贴如下,该类为C++模板类,之所以设计为模板类是为了使得该类支持任意数值类型的卷积计算。其中成员函数ExecuteInputView对应视角一的算法,成员函数ExecuteOutputView对应视角二的算法。想通过代码更加深入理解卷积计算过程的读者可以仔细研究这两个函数的函数体。
#pragma once
#include <vector>
#include <type_traits>
#include <stdio.h>
template <typename T>
class Convolution
{
public:
Convolution();
public:
std::vector<T> ExecuteInputView(const std::vector<T> &x, const std::vector<T> &h);
std::vector<T> ExecuteOutputView(const std::vector<T> &x, const std::vector<T> &h);
private:
T ReverseDotProduct(const T *x, const T *y, int len);
};
template <typename T>
Convolution<T>::Convolution()
{
if (!std::is_arithmetic<T>::value) {
std::exception e("Convolution template class only support arithmetic type.");
throw e;
}
}
template <typename T>
std::vector<T> Convolution<T>::ExecuteInputView(const std::vector<T> &x, const std::vector<T> &h)
{
int nInputSignalNum = static_cast<int>(x.size());
int nPulseResponseNum = static_cast<int>(h.size());
int nOutputSignalNum = nInputSignalNum + nPulseResponseNum - 1;
std::vector<T> y(nOutputSignalNum, static_cast<T>(0));
for (int i = 0; i < nInputSignalNum; ++i) {
for (int j = 0; j < nPulseResponseNum; ++j) {
y[i + j] += x[i] * h[j];
}
}
return y;
}
template <typename T>
std::vector<T> Convolution<T>::ExecuteOutputView(const std::vector<T> &x, const std::vector<T> &h)
{
int nInputSignalNum = static_cast<int>(x.size());
int nPulseResponseNum = static_cast<int>(h.size());
int nOutputSignalNum = nInputSignalNum + nPulseResponseNum - 1;
std::vector<T> xEx(nInputSignalNum + 2 * (nPulseResponseNum - 1), static_cast<T>(0));
memcpy(xEx.data() + nPulseResponseNum - 1, x.data(), sizeof(T)*nInputSignalNum);
std::vector<T> y(nOutputSignalNum, static_cast<T>(0));
for (int i = 0; i < nOutputSignalNum; ++i) {
y[i] = ReverseDotProduct(xEx.data() + i, h.data(), nPulseResponseNum);
}
return y;
}
template <typename T>
T Convolution<T>::ReverseDotProduct(const T *x, const T *y, int len)
{
T res = static_cast<T>(0);
for (int i = 0, j = len - 1; i < len; ++i, --j) {
res += x[i] * y[j];
}
return res;
}
如下代码是卷积计算类的单元测试代码,感兴趣的读者可以在自己的电脑上运行。C++11以上即可使用。
#include <iostream>
#include "Convolution.h"
using namespace std;
int main()
{
vector<int> x = { 1,2,1,2,4,2,1 };
vector<int> h = { 1,2,2,3 };
try {
Convolution<int> conv;
vector<int> yInput = conv.ExecuteInputView(x, h);
vector<int> yOutput = conv.ExecuteOutputView(x, h);
vector<int> ref = { 1, 4, 7, 11, 16, 17, 19, 18, 8, 3 };
if (yInput != ref || yOutput != ref) {
cerr << "UTConvolution failed." << endl;
return 1;
}
}
catch (const exception &e) {
cerr << "UTConvolution exception:" << e.what() << endl;
return 1;
}
cout << "UTConvolution pass." << endl;
1
return 0;
}
以上内容如果存在任何错误和不足,欢迎各位朋友指正,我定会认真接受!!!