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

离散卷积的算法分析及MATLAB实现（实例程序）摘自：张登奇,陈佳.离散卷积的算法分析及MATLAB实现J.湖南理工学院学报(自然科学版). 2013(02)摘 要：离散卷积是信号处理的基本运算,快速卷积和分段卷积是计算离散卷积的重要算法.文章以离散线性卷积的概念为基础,介绍了计算卷积的常用方法和运算流程,列举了MATLAB实现的程序.这些程序能动态演示卷积运算的全部过程,既可帮助理解卷积运算的原理流程,也可作为教学讲解的演示工具.1、离散卷积常规算法的动态演示程序clear;clc;close all;format compact;xn=5,4,3,2,1,M=length(xn),%输入任意序列并计算长度Mhn=1,1,1,N=length(hn),%输入任意脉冲响应并计算长度Nm=-(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,%每算完一位卷积输出后暂停,按任意键继续endyn, %显示线性卷积结果 2、循环卷积运算过程的动态演示程序clc;close all;clear;format compact;xn=5,4,3,2,1,M=length(xn),%输入任意序列并计算长度Mhn=1,1,1,N=length(hn),%输入任意脉冲响应并计算长度NL=8,if Lmax(M,N); L=max(M,N),end,%输入循环卷积的长度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,%每算完一位循环卷积后暂停,按任意键继续endycn, %显示循环卷积结果3、快速算法计算线性卷积的程序clc;close all;clear;format compact;xn=5,4,3,2,1, M=length(xn),%输入任意序列并计算长度Mhn=1,1,1, N=length(hn),%输入任意脉冲响应并计算长度NL=2nextpow2(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),%输入任意序列并计算长度Lxhn=1,1,1,1,N=length(hn),%输入任意脉冲响应并计算长度NL=2nextpow2(2*N),%确定适合基2-FFT算法的长度M=L-N+1,T=ceil(Lx/M),%计算每段的点数M及总段数Lxn=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);%控制绘图坐标ylabel(x(n),grid,title(a) 输入序列x(n),%加标签网格和标题subplot(6,1,2);stem(n,hn,zeros(1,(T+1)*M-N),r.),%绘脉冲响应axis(nmin,nmax,hnmin,hnmax);%控制绘图坐标ylabel(h(n),grid,title(b) 脉冲响应h(n),%加标签网格和标题for i=0:T;disp(分段卷积段号：),i, %依次分段处理并显示段号 ix=i*M+1,%各段在MATLAB变量中的起始位置 xi_seg=xn(ix:ix+M-1),%从输入序列中取出一段 xin=zeros(1,(T+1)*M);%为了只绘第i段输入序列先置零处理 subplot(6,1,3);xin(ix:ix+M-1)=xn(ix:ix+M-1);stem(n,xin,b.), % 绘第i段输入 axis(nmin,nmax,xnmin,xnmax);%控制绘图坐标 ylabel(xi(n),grid,title(c)第i段输入xi(n),%加标签网格和标题 yi_seg=conv(xi_seg,hn),%分段卷积实际上是用快速算法且X(k)在循环前算出保存 yin=zeros(1,(T+1)*M);%为了只绘第i段输出序列先置零处理 if iT; yin(ix:ix+L-1)=yi_seg(1:L);else yin(ix:ix+M-1)=0;end,%最后段输出特殊处理 subplot(6,1,4);stem(n,yin,b.),% 绘第i段输出 axis(nmin,nmax,ynmin,ynmax); %控制绘图坐标 ylabel(yi(n),grid,title(d)第i段输出yi(n),%加标签网格和标题 tin=zeros(1,(T+1)*M);%为了只绘前段缓存序列先置零处理 subplot(6,1,5);tin(ix:ix+N-2)=tn(1:N-1);stem(n,tin,b.),%绘前段缓存序列 axis(nmin,nmax,ynmin,ynmax); %控制绘图坐标 ylabel(t(n),grid,title(e)前段缓存t(n),%加标签网格和标题 yi_seg(1:N-1)=yi_seg(1:N-1)+tn(1:N-1),%重叠相加运算 yn(ix:ix+M-1)=yi_seg(1:M),tn(1:N-1)=yi_seg(M+1:L),%取前M个数据