前言
本人在对野外仪器数据采集的信号做一些列数字信号处理的过程中，经常遇到卷积与循环卷积。现实中，野外仪器记录到的数据都是离散的！所以我们在对数据进行FFT变换时，都默认使用的是DFT(离散傅里叶变换)。
本文目标并非介绍卷积与循环卷积的详细内容与公式，而是侧重如何编写自己的程序实现这两种操作，并简要介绍两者在使用和结果上的区别。
卷积与循环卷积区别
明确说明一点：虽然这两种操作都是对"两个信号"的处理，但两者不是同一操作！
区别1：卷积可以处理连续信号和离散信号，循环卷积只用于处理离散信号(现实中都可用)；
区别2：卷积对a和b这两个信号的长度没有要求，循环卷积要求a和b信号长度一致；
区别3：设操作后的信号为c，卷积操作后c的长度为a+b-1，循环卷积长度为a(a=b)。
卷积如何使用
卷积的操作就是"循环乘积与加和"。下面结合一个具体的例子进行说明：
信号a：a = [1 3 4 9 8 7]
信号b：b = [2 4 6 3]
卷积：a*b
对应的Matlab语句是conv(a,b)。但是用edit conv命令查看conv的源代码时候，感觉写的很复杂，不如根据上面计算的流程图(计算原理)写自己的代码：
a = [1 3 4 9 8 7];
b = [2 4 6 3]; % a*b a在前b在后
a_length = length(a);
b_length = length(b);
r = zeros(1,a_length+b_length-1); % 记录卷积结果，总长度是a+b-1
for k = 1:a_length % a做循环总长 因为a在前所以是a在移动!
c = a(k)*b; % a所有的数都会和b数组所有元素乘一遍，但各自在r中起作用的位置不同
d = r(k:k+b_length-1); % r(k:k+wb-1)很妙！提取当前a(k)元素会影响的r的区间
d = d+c; % 把当前a(k)影响的范围"对应"加进去
r(k:k+b_length-1) = d; % 把r当前受影响的区域更新(加进去)
end
fprintf('卷积结果为:');
r
注意：进行卷积计算的两个信号谁在前谁在后结果都是一样的，即a*b = b*a；但是上面的程序还是要区分一下"谁在前"的问题(小细节，注意一下即可)。
循环卷积如何使用
循环卷积的操作一般用"傅里叶逆变换"来计算。下面结合一个具体的例子进行说明：
信号a：a = [2,1,2,1,6,2,8,9]
信号b：b = [3,4,2,4,3,5,1,8]
循环卷积：a*b
注意：两个信号的长度一样
用自己编的DFT与IDFT进行循环卷积的Matlab程序：
clc; clear;
% 手动输入的信号:
fprintf('原始信号为:\n');
a = [2,1,2,1,6,2,8,9]
b = [3,4,2,4,3,5,1,8]
N = length(a);
WN = exp(-i*2*pi/N); % 常数
WN_nk_a = zeros(N)+WN; % a的WN_kn
WN_nk_b = zeros(N)+WN; % b的WN_kn
xk_a = a'; % a时域信号振幅(列矩阵)
xk_b = b'; % b时域信号振幅(列矩阵)
E_a = zeros(N); % a的辅助的E(WN_kn的幂,单独拿出来算)
E_b = zeros(N); % a的辅助的E(WN_kn的幂,单独拿出来算)
%%% a,b信号傅里叶正变换即结果 %%%
for row = 0:N-1
for cow = 0:N-1
E_a(row+1,cow+1) = row*cow;
E_b(row+1,cow+1) = row*cow;
WN_nk_a(row+1,cow+1) = WN_nk_a(row+1,cow+1)^E_a(row+1,cow+1);
WN_nk_b(row+1,cow+1) = WN_nk_b(row+1,cow+1)^E_b(row+1,cow+1);
end
end
Xk_a = WN_nk_a * xk_a;
Xk_b = WN_nk_b * xk_b;
Xk_t = Xk_a.*Xk_b;
%%% T信号(循环卷积)傅里叶逆变换即结果 %%%
WN_nk_t = zeros(N)+WN; % t的WN_kn
E_t = zeros(N); % t的辅助的E(WN_kn的幂,单独拿出来算)
for row = 0:N-1
for cow = 0:N-1
E_t(row+1,cow+1) = -row*cow;
WN_nk_t(row+1,cow+1) = WN_nk_t(row+1,cow+1)^E_t(row+1,cow+1);
end
end
fprintf('a,b信号的循环卷积结果为:\n');
xk_t = real((WN_nk_t * Xk_t)/N)'
优势1：信号可用是任意长度(不用像FFT那种必须是2^N的长度)；
优势2：完全理解其内涵。
当然，全部用Matlab自带函数写就是下面：
fprintf('原始信号为:\n');
a = [2,1,2,1,6,2,8,9]
b = [3,4,2,4,3,5,1,8]
a_tmp = fft(a);
b_tmp = fft(b); % 用FFT计算a与b信号的频域信号
r_tmp = a_tmp.*b_tmp; % 对应元素相乘(点乘)
r = ifft(r_tmp) % 反变换回去，就是a,b信号循环卷积结果(时域)
