快速卷积的MATLAB实现

摘要
在信号处理中，许多具体的应用是以线性卷积为基础的。

当序列点数较少时可以直接计算线性卷积，然而当序列长度很长时，直接计算卷积的运算量非常庞大。

快速卷积是实现卷积的一种快速算法，减少了运算量，节约了时间，给我们计算卷积提供了很大的便利。

本课程设计是以Matlab为基础，完成序列的卷积和快速卷积运算的编程实现，以及相应的分析和比较。

关键字：Matlab 卷积快速卷积
目录
1.基于设计题目的原理简介 (1)
1.1 序列卷积的定义 (1)
1.2 快速傅里叶变换FFT概念 (1)
1.3 快速卷积方法及实现 (1)
2.程序设计及运行结果分析 (3)
2.1 题目一 (3)
2.2 题目二 (4)
2.3 题目三 (7)
3.心得体会 (10)
参考文献 (11)
专业综合课程设计成绩评定表 (12)
1.基于设计题目的原理简介
卷积是数字信号处理（DSP）系统中最常见的，也是最重要的运算之一，无论在时域或频域都离不开卷积运算，FFT是DFT的快速算法，当满足一定条件时可用来计算线性卷积，称为快速卷积。

Matlab具有强大的矩阵运算能力，方便实用的绘图功能和语言的高度集成性，在DSP开发中，使用Matlab可以快速对系统进行仿真运算。

1.1 序列卷积的定义
设x(n)和h(n)是两个离散序列，进行下列求和运算：
∑∞-∞=
= -
=
n
n
h
n
x
m
n
h
m
x
n
y)
(
*)
(
)
(
)
(
)
(
这样，随着n的不同取值，这个求和公式就定义了一个新序列y(n)，称为序列x(n)与h(n)的卷积，记为y(n)=x(n)*h(n) 。

由于DSP主要依靠计算机完成，而计算机无论在时域或频域只能处理有限长的离散信号。

此时只需令上述公式中的n在一个范围内取值即可。

1.2 快速傅里叶变换FFT概念
DFT就是对序列频谱的离散化，在数字信号处理中有着重要的作用，但直接计算DFT 的运算量非常大，它与序列长度的平方成正比，因此制约了DFT的应用。

快速傅里叶变换FFT是实现DFT的一种快速算法，能使计算机计算离散傅里叶变换所需要的乘法次数大为减少，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著。

因而FFT也有重要的作用，下面一节要介绍的快速卷积就是其应用之一。

1.3 快速卷积方法及实现
在信号处理中，许多具体的应用是以线性卷积为基础的。

我们知道，当满足一定条件(L≥M+N-1)时，可以用圆周卷积来计算线性卷积。

由圆周卷积定理知道，圆周卷积可以借助DFT来运算，因此DFT的快速算法FFT就可以用来计算线性卷积。

设x1(n)与x2(n)分别是长度为N与M的有限长序列，它们的线性卷积为y l(n)，L点的圆周卷积为y c(n)，它们的关系为:
y c(n)=∑y l(n+rL)R L(n)
快速卷积算法过程的示意图如图1.1所示。

在实际应用中，常遇到的问题是参加卷积的两个序列的长度相差较大，这样长度小的序列就需补很多的零点，这样就需要大的存储量，运算时间也会变长。

常用的解决方法有两种，一是重叠想加法，另一种是重叠保留法。

这里不作介绍。

图1.1 快速卷积示意图
1.4 直接卷积和快速卷积分析比较
快速卷积，顾名思义，其重点在一个“快”，如果对卷积速度要求较高，快速卷积无疑是理想的工具。

当然，为了提高速度，就要牺牲面积和功耗。

而且由上一节的介绍可知，快速卷积运算的步骤较多，因而当序列较短时快速卷积运算并没有优势，直接进行卷积运算反而更简便。

卷积计算的方法选择要视实际情况而定。

2.程序设计及运行结果分析
2.1 题目一
已知线性非移变系统的h(n)=[6,2,3,6,4,2]，输入为x(n)=[1,2,3,4,5]；
（1）用人工计算系统输出y(n)；
（2）编写程序输出y(n)，并作图。

2.1.1 人工计算
长度为N=6的序列h(n)和长度为M=5的序列x(n)，卷积y(n)的序列长度为(M+N-1)=10，计算过程如下：
y(1) = h(1)•x(1) =6
y(2) = h(1)•x(2)+h(2)•x(1) =14
y(3) =h(1)•x(3)+h(2)•x(2)+h(3)•x(1) = 25
…
y(n) = h(1)•x(n)+h(2)•x(n-1)+ … +h(n)•x(1) …
y(M+N-1) = h(1)•x(M+N-1)+h(2)•x(M+N)+h(3)•x(M+N+1)+…+h(M+N-1)•x(1) =10 所得的结果为y(n)= [6，14，25，36，63，50，55，52，28，10]。

2.1.2 程序设计
在Matlab中实现卷积的函数conv，可以直接调用。

设计程序如下：
h=[6,2,3,6,4,2];
x=[1,2,3,4,5];
y=conv(h,x); %调用conv函数直接计算线性卷积
stem(y); %画出卷积结果h(n)的序列图
title('y(n)')
2.1.3 运行结果及分析
程序运行后，得到的卷积结果y(n)如图2.1所示。

由图2.1
设计的程序是正确的，得到了正确的结果。

图2.1 编程得到的卷积结果y(n)
2.2 题目二
用函数conv和FFT计算长为1000序列的卷积，比较其计算时间。

2.2.1 设计内容及原理分析
本题目的在于比较直接卷积和快速卷积的计算时间，分析其优劣性。

由1.3节可知快速卷积算法如下：
（1）序列补零：将两序列都补零到L点；
（2）计算X1(k)=FFT[x1(n)]；
（3）计算X2(k)=FFT[x2(n)]；
（4）计算Y(k)=X1(k)X2(k)；
（5）计算y(n)=x1(n)*x2(n)=IFFT[Y(k)]。

本题中已经给定两原序列的长度M=N=1000，则L≥M+N-1，而因为快速卷积的基础是FFT，所以要求L满足L=2n，为减少运算量选取L=2048。

设两原序列分别为：
x1=0.5sin(2*n) ；x2=n^3
计算时间的获取方法为：计算前先调用clock函数读取瞬时时钟，待计算结束后，调用etime(t1，t2)函数计算时刻t1，t2间所经历的时间。

2.2.2 程序设计框图
图2.2 程序设计框图
2.2.3 程序代码
N=1000;
L=pow2(nextpow2(1000+1000-1));%计算L的值
n=1:N;
x1=0.5*sin(2*n); x2=n.^3;
t0=clock; %调用clock函数读取瞬时时钟
yc=conv(x1,x2); %用函数conv计算卷积
conv_time=etime(clock,t0) %调用etime函数计算时间
t0=clock;
yf=ifft(fft(x1,L).*fft(x2,L)); %用函数FFT计算卷积
fft_time=etime(clock,t0)
subplot(321),stem(x1,'.');ylabel('x1(n)');
subplot(322),stem(x2,'.');ylabel('x2(n)');
subplot(312);stem(real(yc),'.');ylabel('直接卷积y(n)');
subplot(313);stem(real(yf),'.');ylabel('快速卷积y(n)');
2.2.4 运行结果及分析
用函数conv和FFT计算长为1000序列的卷积的结果如图2.3，计算时间如图2.4。

可以看出，两种方法所计算出的卷积结果是一样的。

再来看时间，用函数conv计算该卷积需要0.0070s，而用FFT计算该卷积需要0s。

由此可以得出，当序列长度为1000时，FFT算法所消耗的时间远小于函数conv要消耗的时间。

图2.3 函数conv和FFT计算长为1000序列的卷积的结果
图2.4 函数conv 和FFT 的计算时间
2.3 题目三
用快速卷积法计算()0.9()n M x n R n =和()()N h n R n =两个序列的卷积；并测试直接卷积和快速卷积的时间。

2.3.1 设计内容及原理分析
用快速卷积计算两个序列卷积的方法与题目二中的一样，区别仅在于序列的长度没有给出，这里不再重复计算过程。

为了做对比分析，两个原序列的长度M 和N 取两组值进行仿真，分别取M=50，N=46和M=800，N=1000。

2.3.2 程序设计框图
程序设计框图与题目二一样，如图2.2所示。

2.3.3 程序代码
xn=0.9.^(1:800);
hn=1.^(1:1000);
L=pow2(nextpow2(800+1000-1)); %计算L 值
tic;
yc=conv(xn,hn); %直接计算卷积
toc; %计算时间
tic;
Xk=fft(xn,L); %求x(n)的快速傅里叶变换X(k)
Hk=fft(hn,L); %求h(n)的快速傅里叶变换H(k)
Yk=Xk.*Hk; %求Y(k)
yf=ifft(Yk,L); %用IFFT求快速卷积yf
toc; %计算时间
subplot(221),stem(xn,'.');title('x(n)');
subplot(222),stem(hn,'.');title('h(n)');
subplot(212),ny=1:L;stem(real(yf),'.');title('快速卷积y(n)');
2.3.4 运行结果及分析
（1）M=50，N=46
卷积结果及计算时间分别如图2.5和图2.6所示（直接卷积的计算时间在上，快速卷积的计算时间在下）。

图2.5 M=50，N=46时的快速卷积
图2.6 M=50，N=46时的计算时间
（2）M=800，N=1000
卷积结果及计算时间分别如图2.7和图2.8所示（直接卷积的计算时间在上，快速卷积的计算时间在下）。

图2.7 M=800，N=1000时的快速卷积
图2.8 M=800，N=1000时的计算时间
可以看出，当M=50，N=46时，快速卷积的计算时间比直接卷积的要长；而当M=800，N=1000时，快速卷积的时间更短。

这说明：当序列的点数比较少时，快速卷积并不占有优势，所用时间反而比直接卷积长，只有在序列较长时，快速卷积才体现出“快”的优势。

3.心得体会
本次课程设计我的题目是——序列的卷积和快速卷积运算的编程实现，因为我们已经学习过信号与系统和Matlab应用实践课程，对序列的卷积和Matlab仿真并不陌生；而且此次题目卷积和快速卷积在课本上能找到相关的内容，老师在授课时已经讲得很详细，我对这部分的内容也比较熟悉。

所以这次课程设计做的比较顺利，没有遇到什么大的难题，只是有一些小的细节影响了设计的进程，比如函数名的拼写等。

这让我体会到，即使是比较容易的事情时，在对待它时也要很认真仔细，就像我们常说的：细节决定成败。

数字信号处理是本专业的一门重要课程，与它相关性最强的软件就是Matlab，这次课程设计给我提供了一个用Matlab实现序列卷积和快速卷积的机会，通过编程得到了和理论一致的结果。

不仅加深了我对序列卷积和快速卷积这部分内容的理解，而且提高了我的程序编写和纠错能力，也进一步熟悉了Matlab的操作。

但这是远远不够的，卷积只是数字信号处理中比较简单的一部分，这门课程还有许多值得我们去实践和探究的地方，我现在的目标就是在假期将本课程其它几个课设题目都操作一遍，以进一步提升自己各方面的能力。

作为一名电子信息工程专业的学生，熟练掌握相关软件的操作是一门必修课，否则就算理论知识学得再好，不会应用，也只是纸上谈兵，不会真正取得什么成果，毕业之后在工作上也会很吃力。

而我们并没有专门的课程来学习这些软件，所以这就需要我们平时多多利用课余时间来自学。

现在互联网又发达，我们应充分利用这些资源，通过不断地积累知识来提高的能力。

参考文献
[1]刘泉等.数字信号处理原理与实现（第2版）[M].北京：电子工业出版社，2009.
[2] 程佩青.数字信号处理教程[M].北京：清华大学出版社，2008.
[3]李建新等.现代通信系统分析与仿真-MATLAB通信工具箱[M].西安：西安电子科技大学出版社，2000.
[4]丁玉美等.数字信号处理[M].西安：西安电子科技大学出版社，2001.
[5]陈怀琛等.MATLAB及在电子信息课程中的应用[M].北京：电子工业出版社，2003.。