信号处理--线性卷积与循环卷积

一、线性卷积

线性卷积（linear convolution）是在时域描述线性系统输入和输出之间关系的一种运算。

1. 应用背景

  对于线性时不变离散时间系统来说，若序列$x(n)$是系统的输入，$h(n)$是系统在单位脉冲作用下的单位脉冲响应，由于输入离散时间序列$x(n)$可表示为一系列脉冲的线性组合，根据线性系统的齐次性和可加性，$x(n)$作用于系统所引起的零状态响应$y(n)$就是序列$x(n)$与$h(n)$的卷积和。

2. 定义式

$$y(n)=x(n)\ast h(n)=\sum _{k=-\infty }^{\infty }x(k)h(n-k)$$

上式在运算过程中存在序列的翻转、移位、相乘和相加，所以称为卷积和。

3. 计算方法

3.1 定义式

3.2 作图法

• 步骤：

（1）将序列$x(n)$、$h(n)$的自变量用$k$代换，然后将序列$h(k)$以纵坐标为轴线反转，成为$h(-k)$；
（2）序列$h(-k)$沿$k$轴正方向平移$n$个单位，成为$h(n-k)$；
（3）求乘积$x(k)h(n-k)$；
（4）按定义式求各乘积之和。

3.3 列表法

由$$y(n)=x(n)\ast h(n)=\sum _{k=-\infty }^{\infty }x(k)h(n-k)$$可见，求和符号内$x(k)$的序号$k$与$h(n-k)$的序号$n-k$之和`恰好等于$k$。

沿斜线上各数值之和就是卷积和。

---

二、循环卷积

1. 序列的循环移位

设$x(n)$为有限长序列，长度为M，M≤N，则$x(n)$的循环移位定义为：
$$y(n)=x((n+m){)}_N{R}_N(n)$$
该式表明，将$x(n)$以N为周期进行周期延拓得到$\tilde{x}(n)=x((n){)}_N$，再将$\tilde{x}(n)$左移m得到$\tilde{x}(n+m)=x((n+m){)}_N$，最后取$\tilde{x}(n+m)$的主值序列则得到有限长序列$x(n)$的循环移位序列$y(n)$。$y(n)$是长度为N的有限长序列。

clear all;
close all;
clc;
%%%%%%%%%%%%%%% x(n)及其循环移位过程 %%%%%%%%%%%%%
M=24;
xn=triang(6);
n=-M/3:2*M/3-1;
x=zeros(1,M);
for l=1:6
    x(8+l)=xn(l);
end
subplot(4,2,1);
stem(n,x,'.'); axis([-8,16,0,1]); title('(a) x(n)');
b=[xn' 0 0];
y8=[b,b,b];
subplot(4,2,3);
stem(n,y8,'.'); axis([-8,16,0,1]); title('(b) x((n))_8')
y8=[y8 b];
for l=1:M'
    ym(l)=y8(l+2);
end
subplot(4,2,5);
stem(n,ym,'.'); axis([-8,16,0,1]); title('(c) x((n+2))_8');
for l=9:16
    x(l)=ym(l);
end
subplot(4,2,7);
stem(n,x,'.'); axis([-8,16,0,1]); xlabel('n'); title('(d) x((n+2))_8R_8(n)');

2. 循环卷积的定义

设序列$h(n)$和$x(n)$的长度分别为N和M，$h(n)$与$x(n)$的L点循环卷积定义为：
$$y_c(n)=[\sum _{m=0}^{L-1}h(m)x((n-m){)}_L]R_L(n)$$
式中，L称为循环卷积区间长度，$L\ge max[N,M]$，$x((n-m){)}_L$是以L为周期的周期信号，n和m的变化区间均是$[0，L-1]$。
因此直接计算该式比较麻烦，计算机中采用矩阵相乘或FFT的方法计算循环卷积。

3. 用矩阵计算循环卷积的公式

  当$n=0,1,2,\cdots ,L-1$时，由$x(n)$形成的序列为 { x ( 0 ) , x ( 1 ) , x ( 2 ) , ⋯   , x ( L − 1 ) } \{x(0),x(1),x(2),\cdots,x(L-1)\} {x(0),x(1),x(2),⋯,x(L−1)}。

  令$n=0,m=0,1,2,\cdots ,L-1$，$x((n-m){)}_L$形成的序列为：
{ x ( ( 0 ) ) L , x ( ( − 1 ) ) L , x ( ( − 2 ) ) L , ⋯   , x ( ( − L + 1 ) ) L } = { x ( 0 ) , x ( L − 1 ) , x ( L − 2 ) , ⋯   , x ( 1 ) } \begin{aligned} &\{x((0))_L,x((-1))_L,x((-2))_L,\cdots,x((-L+1))_L\}\\ &=\{x(0),x(L-1),x(L-2),\cdots,x(1)\} \end{aligned} ​{x((0))L​,x((−1))L​,x((−2))L​,⋯,x((−L+1))L​}={x(0),x(L−1),x(L−2),⋯,x(1)}​
与序列$x(n)$进行对比，相当于将第一个序列值$x(0)$保持不变，将后面的序列 { x ( 1 ) , x ( 2 ) , ⋯   , x ( L − 1 ) } \{x(1),x(2),\cdots,x(L-1)\} {x(1),x(2),⋯,x(L−1)}反转180°再放到$x(0)$的后面。这样形成的序列称为$x(n)$的循环倒相序列。

  令$n=1,m=0,1,2,\cdots ,L-1$，$x((n-m){)}_L$形成的序列为：
{ x ( ( 1 ) ) L , x ( ( 0 ) ) L , x ( ( − 1 ) ) L , ⋯   , x ( ( − L + 2 ) ) L } = { x ( 1 ) , x ( 0 ) , x ( L − 1 ) , ⋯   , x ( 2 ) } \begin{aligned} &\{x((1))_L,x((0))_L,x((-1))_L,\cdots,x((-L+2))_L\}\\ &=\{x(1),x(0),x(L-1),\cdots,x(2)\} \end{aligned} ​{x((1))L​,x((0))L​,x((−1))L​,⋯,x((−L+2))L​}={x(1),x(0),x(L−1),⋯,x(2)}​
观察上式等号右端序列，它相当于$x(n)$的循环倒相序列向右循环移一位，即向右移1位，移出区间$[0，L-1]$的序列值再从左边移进。

依次类推，当 n 和 m 均从 0 变化到 L-1 时，得到一个关于$x((n-m){)}_L$的矩阵：
$$\left[\begin{array}{ccccc}x(0)& x(L-1)& x(L-2)& \cdots & x(1)\\ x(1)& x(0)& x(L-1)& \cdots & x(2)\\ x(2)& x(1)& x(0)& \cdots & x(3)\\ ⋮& ⋮& ⋮& ⋮& ⋮\\ x(L-1)& x(L-2)& x(L-3)& \cdots & x(0)\end{array}\right]$$
上面矩阵称为x(n)的L点“循环卷积矩阵”，其特点是：

1. 第一行是序列 { x ( 0 ) , x ( 1 ) , x ( 2 ) , ⋯   , x ( L − 1 ) } \{x(0),x(1),x(2),\cdots,x(L-1)\} {x(0),x(1),x(2),⋯,x(L−1)}的循环倒相序列。注意，如果$x(n)$的长度M<L，则需要在$x(n)$末尾补 L-M 个零后，再形成第一行的循环倒相序列；
2. 第一行以后的各行均是前一行向右循环移1位形成的；
3. 矩阵的各主对角线上的序列值均相等。

有了循环卷积矩阵，就可以写出计算循环卷积的矩阵形式：

$\left[\begin{array}{c}{y}_c(0)\\ y_c(1)\\ y_c(2)\\ ⋮\\ y_c(L-1)\end{array}\right]=\left[\begin{array}{ccccc}x(0)& x(L-1)& x(L-2)& \cdots & x(1)\\ x(1)& x(0)& x(L-1)& \cdots & x(2)\\ x(2)& x(1)& x(0)& \cdots & x(3)\\ ⋮& ⋮& ⋮& ⋮& ⋮\\ x(L-1)& x(L-2)& x(L-3)& \cdots & x(0)\end{array}\right]\left[\begin{array}{c}h(0)\\ h(1)\\ h(2)\\ ⋮\\ h(L-1)\end{array}\right]$

上式中如果$h(n)$的长度N<L，则需要在$h(n)$末尾补上 L-N 个零。

三、线性卷积和循环卷积之间的关系

• 假设$h(n)$和$x(n)$都是有限长序列，长度分别为N和M。它们的线性卷积和循环卷积分别表示为：
  $$y_l(n)=h(n)\ast x(n)=\sum _{m=0}^{N-1}h(m)x(n-m)$$
  $$y_c(n)=[\sum _{m=0}^{L-1}h(m)x((n-m){)}_L]R_L(n)$$
  其中，$L\ge max[N,M]$，$x((n){)}_L={\sum }_{i=-\infty }^{\infty }x(n+iL)$
  所以$$\begin{array}{rl}{y}_c(n)& =[\sum _{m=0}^{N-1}h(m)\sum _{i=-\infty }^{\infty }x(n-m+iL)]R_L(n)\\ & =[\sum _{i=-\infty }^{\infty }\sum _{m=0}^{N-1}h(m)x(n+iL-m)]R_L(n)\end{array}$$
  与$y_l(n)$式对照可以看出
  $$\sum _{m=0}^{N-1}h(m)x(n+iL-m)=y_l(n+iL)$$
  即$$y_c(n)=[\sum _{i=-\infty }^{\infty }{y}_l(n+iL)]R_L(n)$$
  也就是说，$y_c(n)$等于$y_l(n)$以L为周期的周期延拓序列的主值序列。

• 线性卷积与循环卷积相等的条件

• $y_l(n)$的长度为N+M-1,因此只有当循环卷积长度L≥N+M-1时，$y_l(n)$以L为周期进行周期延拓时才无时域混叠现象，此时取其主值序列显然满足 $y_c(n)=y_l(n)$。

• 由此说明了循环卷积等于线性卷积的条件是L≥N+M-1。

clear all;
close all;
clc;
%%%%%%%%%%%%%%% 验证线性卷积与循环卷积直接的关系 %%%%%%%%%%%%%
xn=ones(1,5);                 % 设置x(n)=R5(n)
hn=ones(1,4);                 % 设置h(n)=R4(n)
Lx=length(xn);
Lh=length(hn);

nh=0:Lh-1;                  
subplot(3,2,1); 
stem(nh,hn,'.'); axis([0,10,0,1.5]);ylabel('h(n)');xlabel('n');title('(a) h(n)波形');

nx=0:Lx-1;
subplot(3,2,3);
stem(nx,xn,'.');axis([0,10,0,1.5]);ylabel('x(n)');xlabel('n');title('(b) x(n)波形');
% 调用conv计算线性卷积
yn=conv(hn,xn);
L=length(yn);
ny=0:L-1;
subplot(3,2,5);
stem(ny,yn,'.');axis([0,10,0,5]);ylabel('y_l(n)');xlabel('n');title('(c) y_l(n)=h(n)*x(n)');
% 循环卷积等于线性卷积的条件:循环卷积区间长度L>=M+N-1
% L=6,用DFT计算序列循环卷积
L1=6;
Xk=fft(xn,L1);                    % L1点FFT[x(n)]
Hk=fft(hn,L1);                    % L1点FFT[h(n)]
Yck=Xk.*Hk;                       % 频域相乘得Yc(k)
ycln=ifft(Yck,L1);                % L1点IFFT得到6点循环卷积结果yc(n)
nycl=0:L1-1;
subplot(3,2,2);
stem(nycl,ycln,'.');axis([0,10,0,5]);ylabel('y_c(n)');xlabel('n');title('(d) 6点循环卷积');
%L=8,用DFT计算序列循环卷积=线性卷积
L2=8;
Xk=fft(xn,L2);                      % L2点FFT[x(n)]
Hk=fft(hn,L2);                      % L2点FFT[h(n)]
Yck=Xk.*Hk;                         % 频域相乘得Yc(k)
yc2n=ifft(Yck,L2);                  % L2点IFFT得到8点循环卷积结果yc(n)
nyc2=0:L2-1;
subplot(3,2,4);
stem(nyc2,yc2n,'.');axis([0,10,0,5]);ylabel('y_c(n)');xlabel('n');title('(e) 8点循环卷积');
%L=10,用DFT计算序列线性卷积
L3=10;
Xk=fft(xn,L3);                      % L3点FFT[x(n)]
Hk=fft(hn,L3);                      % L3点FFT[h(n)]
Yck=Xk.*Hk;                         % 频域相乘得Yc(k)
yc3n=ifft(Yck,L3);                  % L2点IFFT得到10点循环卷积结果yc(n)
nyc3=0:L3-1;
subplot(3,2,6);
stem(nyc3,yc3n,'.');axis([0,10,0,5]);ylabel('y_c(n)');xlabel('n');title('(f) 10点循环卷积');

四、用DFT计算循环卷积

设$h(n)$和$x(n)$的长度分别为N和M，其L点循环卷积为
$$y_c(n)=[\sum _{m=0}^{L-1}h(m)x((n-m){)}_L]R_L(n)$$
且
$$\{\begin{array}{l}H(k)=DFT[h(n){]}_L\\ X(k)=DFT[x(n){]}_L\end{array}0\le k\le L-1，L\ge max[N,M]$$
则由DFT的时域循环卷积定理有
$$Y_c(k)=DFT[y_c(n){]}_L=H(k)X(k)0\le k\le L-1$$
由此可见，循环卷积既可以在时域直接计算，也可以按照下面的计算框图在频域计算。
由于DFT有快速算法，当L很大时，在频域计算循环卷积的速度快得多，因而常用DFT（FFT）计算循环卷积。

DFT只能直接用来计算循环卷积，然而在实际应用中，为了分析时域离散线性时不变系统或者对序列进行滤波处理等，需要计算两个序列的线性卷积。与计算循环卷积一样，为了提高运算速度，也希望用DFT（FFT）计算线性卷积。

五、用DFT计算线性卷积

取L≥N+M-1，则可按照上述框图用DFT（FFT）计算线性卷积。

六、MATLAB实现

1. 利用MATLAB的conv()函数计算线性卷积

% x1的长度N=5，x2的长度M=3,那么x1与x2线性卷积的输出序列长度L=N+M-1=7
% 可以验证，输出结果与我们使用列表法得到的结果一致
x1=[1,1,3,4,2];
x2=[1,3,2];
y=conv(x1,x2);

2. 利用矩阵相乘原理计算循环卷积

function [y] = circle_convolution(x1,x2,L)
% circle_convolution 函数用于计算两个有限长序列的L点循环卷积
% x1: 长度为N的有限长序列
% x2: 长度为M的有限长序列
% L:  循环卷积的长度
% y:  x1与x2的L点循环卷积结果，L≥max(N,M);

%********** 构造L点循环卷积矩阵 **********%
N=length(x1);
M=length(x2);

x_1=[x1 ,zeros(1,L-N)] ;
x_2=[x2 ,zeros(1,L-M)] ;
x=zeros(L,L);
h=circshift ( fliplr (x_1) ,1);  % 得到x1(n)的循环倒相序列
for i=0:L-1
    x(i+1,:)= circshift(h,i);    % x即为L点循环卷积矩阵
end
%********** 矩阵相乘计算循环卷积 **********%
x_2=x_2';     % 将补零后的x2转置为列向量
y=x*x_2;      % 矩阵相乘
end

% 调用circle_convolution函数计算循环卷积
x1=[1,1,3,4,2];
x2=[1,3,2];
L=7;
y=circle_convolution(x1,x2,L);

当L=9时，结果为：

当L=5时，结果为：

这里也进一步验证了线性卷积和循环卷积的关系：
(1) 当L<N+M-1时，循环卷积是线性卷积长度为L的混叠；
(2) 当L=N+M-1时，循环卷积=线性卷积；
(3) 当L>N+M-1时，循环卷积是线性卷积末尾补L-(N+M-1)个零；

3. 调用fft()函数计算循环卷积与线性卷积

% 5点DFT
x1=[1,1,3,4,2];
x2=[1,3,2];

X1=fft(x1,5);
X2=fft(x2,5);
Y=X1.*X2;
y=ifft(Y);

% 7点DFT
x1=[1,1,3,4,2];
x2=[1,3,2];

X1=fft(x1,7);
X2=fft(x2,7);
Y=X1.*X2;
y=ifft(Y);

% 9点DFT
x1=[1,1,3,4,2];
x2=[1,3,2];

X1=fft(x1,9);
X2=fft(x2,9);
Y=X1.*X2;
y=ifft(Y);