离散卷积的算法分析及MATLAB实现(程序实例)

离散卷积的算法分析及MATLAB实现（实例程序）

摘自：张登奇,陈佳.离散卷积的算法分析及MATLAB实现[J].湖南理工学院学报(自然科学版). 2013(02)

摘要：离散卷积是信号处理的基本运算,快速卷积和分段卷积是计算离散卷积的重要算法.文章以离散线性卷积的概念为基础,介绍了计算卷积的常用方法和运算流程,列举了MATLAB实现的程序.这些程序能动态演示卷积运算的全部过程,既可帮助理解卷积运算的原理流程,也可作为教学讲解的演示工具. 1、离散卷积常规算法的动态演示程序

clear;clc;close all;format compact;

xn=[5,4,3,2,1],M=length(xn),%输入任意序列并计算长度M

hn=[1,1,1],N=length(hn),%输入任意脉冲响应并计算长度N

m=[-(M-1):M+N-2],%设置代换变量的范围以便x(m)翻转和移位

xm=[zeros(1,M-1),xn,zeros(1,N-1)],%补零以便与m对应绘图

subplot(4,1,1);stem(m,xm,'r.'),%%绘输入序列x(m)

ylabel('x(m)'), grid,title('(a) 输入序列x(m)'),%%加标签网格和标题

hm=[zeros(1,M-1),hn,zeros(1,M-1)],%补零以便与m对应绘图

subplot(4,1,2);stem(m,hm,'r.'),%%绘脉冲响应

ylabel('h(m)'),grid,title('(b) 脉冲响应h(m)'),%%加标签网格和标题

yn=zeros(1,2*M+N-2);%卷积输出初始化

for n=0:M+N-2,%逐个计算卷积输出

if n==0;xmfy=[fliplr(xn),zeros(1,M+N-2)],%实现翻转

else for k=M:-1:1; xmfy(k+n)=xmfy(k+n-1);end ; xmfy(n)=0; xmfy,%实现翻转后移位并显示

end

subplot(4,1,3);stem(m,xmfy,'b.'),%%绘制翻转移位序列

ylabel('x(n-m)'), grid,title('(c) 翻转移位x(n-m)'),%%加标签网格和标题

yn(M+n)=sum(xmfy.*hm);%计算第n位输出并与m位置对应

subplot(4,1,4);stem(m,yn,'r.'),%%绘制卷积输出序列

axis([min(m),max(m),min([0,conv(xn,hn)]),max([0,conv(xn,hn)])]),%%控制绘图坐标

ylabel('y(n)'), grid,title('(d) 卷积输出y(n)'),%%加标签网格和标题

pause,%每算完一位卷积输出后暂停,按任意键继续

end

yn, %显示线性卷积结果

2、循环卷积运算过程的动态演示程序

clc;close all;clear;format compact;

xn=[5,4,3,2,1],M=length(xn),%输入任意序列并计算长度M

hn=[1,1,1],N=length(hn),%输入任意脉冲响应并计算长度N

L=8,if L

m=0:L-1,%设置代换变量并与两圆周上的编号对应

xm=[xn,zeros(1,L-M)],%按圆周编号方向排列不足补零

subplot(4,1,1);stem(m,xm,'r.'),%%绘输入序列x(m)

ylabel('x(m)'),grid,title('(a) 输入序列x(m)'),%%加标签网格和标题

hm=[hn,zeros(1,L-N)],%按圆周编号方向排列不足补零

subplot(4,1,2);stem(m,hm,'r.'),%%绘脉冲响应h(m)

ylabel('h(m)'),grid,title('(b) 脉冲响应h(m)'),%%加标签网格和标题

ycn=zeros(1,L);%循环卷积输出初始化

for n=0:L-1%逐个计算循环卷积

xn_m=xm(mod(n-m,L)+1),%按圆周编号方向反排后正转n位

subplot(4,1,3);stem(m,xn_m,'b.'),%%绘反排正转序列

ylabel('x((n-m))L'),grid,title('(c) 反排正转x((n-m))L'),%%加标签网格和标题

ycn(n+1)=sum(xn_m.*hm);%计算第n位循环卷积输出

subplot(4,1,4);stem(m,ycn,'r.'),%%绘循环卷积输出序列

axis([min(m),max(m),min([0,conv(xn,hn)]),max([0,conv(xn,hn)])]) %%控制绘图坐标

ylabel('yc(n)'),grid,title('(d) 循环卷积yc(n)'),%%加标签网格和标题

pause,%每算完一位循环卷积后暂停,按任意键继续

end

ycn, %显示循环卷积结果

3、快速算法计算线性卷积的程序

clc;close all;clear;format compact;

xn=[5,4,3,2,1], M=length(xn),%输入任意序列并计算长度M

hn=[1,1,1], N=length(hn),%输入任意脉冲响应并计算长度N

L=2^nextpow2(M+N-1),%确定适合基2-FFT算法的长度

xnL=[xn,zeros(1,L-M)], hnL=[hn,zeros(1,L-N)],%序列补零至L长

Xk=fft(xnL); Hk=fft(hnL);%快速傅里叶变换

Yk=Xk.*Hk; ycn=ifft(Yk),%对应相乘后求逆输出

yn=conv(xn,hn),%与直接线性卷积结果对比

4、重叠相加法进行线性卷积的动态演示程序

clc;close all;clear;format compact;

xn=[4,3,2,1,1,2,3,4,4,3,2,1,1,2],Lx=length(xn),%输入任意序列并计算长度Lx

hn=[1,1,1,1],N=length(hn),%输入任意脉冲响应并计算长度N

L=2^nextpow2(2*N),%确定适合基2-FFT算法的长度

M=L-N+1,T=ceil(Lx/M),%计算每段的点数M及总段数L

xn=[xn,zeros(1,(T+1)*M-Lx)],%输入补零至最后缓存数据可以输出

yn=zeros(1,(T+1)*M);%卷积输出变量初始化

tn=zeros(1,N-1);%缓存变量初始化

n=0:(T+1)*M-1,%序列号赋值

nmin=min(n)-0.01;nmax=max(n)+0.01;%%确定绘图横坐标并保证能绘出首尾数据

xnmin=min([0,xn]);xnmax=max([0,xn]);%%确定输入序列绘图纵坐标

hnmin=min([0,hn]);hnmax=max([0,hn]);%%确定脉冲响应绘图纵坐标

ynmin=min([0,conv(xn,hn)]);ynmax=max([0,conv(xn,hn)]); %%确定输出序列绘图纵坐标subplot(6,1,1);stem(n,xn,'r.'),%%绘输入序列x(n)

axis([nmin,nmax,xnmin,xnmax]);%%控制绘图坐标