基于MATLAB离散卷积的实现和源代码

2离散卷积的基本原理
2.1线性卷积
首先了解到两个连续时间信号x（t）和h（t）的卷积运算为
并在此基础上得道了LTI连续系统零状态响应的时域计算公式。类似的我们定义：
为序列x(n)和h(n)卷积和运算，即就是离散的线性卷积。
2.2 循环卷积
设序列h(n)和x(n)的长度分别为N和M。h(n)与x(n)的L点循环卷积定义为：
signal1=5*sin(2*pi*f*n)+10*cos(2*pi*n/8);
figure(1);
subplot(411)
stem(n,signal1);
title('原信号1');
xlabel('n');ylabel ('y(n)');
axis([0 18 -25 25])
long_M=5;
signal2=2*ones(1,long_M);
在序列比较长时FFT是一种最合适的方法，可以快速得出仿真运算的结果。选择圆周卷积或FFT运算线性离散卷积可提高DSP设计的效率。
通过对离散卷积的实现，加深了对卷积的理解，卷积的实现方式，卷积和FFT对应的关系。
附件:
clc;
clear;
close all;
N=14;
n=[1:N-1];
f=1/16;
if (i>long_M) %判断翻转序列是否超过了自身的长度
fk(k)=a; %超过了，后面的数为零，对应点上的值保持不变
else
fk(k)=a+signal2(i)*signal1(k-i+1);%计算对应点上的卷积值
a=fk(k);
end
end
else
for i=1:1:k %翻转列超过了原序列右边的第最后一个数
grid on;
subplot(412)
stem(circy);
title('用函数实现的循环卷积条件满足L>=M+N-1');
xlabel('n');ylabel('幅度');
grid on;
axis([0 18 -100 100]);
%%%%%%%循环卷积-当长度条件满足L<M+N-1（和线性卷积不相等）%%%
fn=1:1:changdu;
subplot(414)
stem(fn,Cy);
title('用fft方法实现的卷积');
xlabel('(n)');ylabel('幅度');
axis([0 18 -100 100]);
grid on;
xlabel('n');ylabel('幅度');
axis([0 18 -100 100]);
grid on;
%%%%%%%%%%fft实现卷积%%%%%%%%%%%%%
X1=fft(signal1,changdu);
X2=fft(signal2,changdu);
Y=X1.*X2;
Cy=ifft(Y,changdu);
xlabel('n');ylabel('幅度');
axis([0 18 -100 100]);
grid on;
%%%%%%%%循环卷积-当长度条件满足L>=M+N-1%%%%%%%%%%%%%%%
changdu=long_N+long_M-1;
Blsignal1=[signal1,zeros(1,changdu-long_N)];%补零
5分析的结论和心得体会
通过对线性卷积函数的重新实现后得到的序列，验证了函数是正确的，但是由于采用了较多for循环，使得此函数的运行效率不是很高。
用循环卷积和FFT运算离散卷积具有很大优势，在Matlab中循环卷积没有专用函数，需要自己编制，此例分析了循环卷积条件对于结果的影响，第三个序列是条件不满足N>=N1+N2-1时，序列变短且值对应的值也不同。
Blsignal2=Hn(L1,:); %移位赋值
end
end
circy=Hn'*Blsignal1';
figure(2);
subplot(411)
stem(JuanJi);
title('调用conv函数实现的卷积');
xlabel('n');ylabel('幅度');
axis([0 18 -100 100])
fk=zeros(0,loபைடு நூலகம்g_N+long_M+10); %初始化一个矩阵fk
if(long_N>long_M) %选择较短的矩阵进行翻转的运算
for k=1:1:long_N+long_M-1
a=0;
if (k<=long_N)%判断翻转序列的右边第一个数是否超过了原序列的最后一个数
for i=1:1:k
] (n)
L为循环卷积区间长度，L>=max[N,M]。本为的程序实现方法采用的是矩阵相乘实现循环卷积，
2.3 用FFT实现卷积
根据信号在时域与频域的关系，可以用FFT计算线性卷积，有定义可知
可以对x(n)和h(n)做N=N1+N2-1的FFT（其中N1为x(n)的长度，N2为h(n)的长度），再对 做N点得IFFT即可得到y(n)。
subplot(413)
stem(JuanJi);
title('用conv函数实现的卷积');
xlabel('n');ylabel('幅度');
axis([0 18 -100 100]);
grid on;
%%%%%%%%%%%重新对函数的实现%%%%%%%%%
long_N=length(signal1);
3离散卷积分析的步骤
1.通过定义重新实现线性卷积函数，对比调用的conv（）函数验证重新实现卷积的是否正确。
2.实现循环卷积，通过设定循环卷积的长度条件N>=N1+N2-1和N<N1+N2-1，分析循环卷积和线性相等时满足的条件。
3.将两个序列做N点得FFT，举证相乘后，再做N点得IFFT，从而得到卷积的序列。
4离散卷积分析的结果
图1
图1为Matlab的程序仿真图，第三个序列为调用函数conv()实现的线性卷积，第四个序列为自编线性卷积函数。可以从两图中看到序列是相同的，从而说明自编的线性卷积函数是正确的。
图2
图2中的第一个是线性卷积的实现，第二个序列是条件满足N>=N1+N2-1时的循环卷积，第三个序列是条件不满足N>=N1+N2-1时的循环卷积，第四个序列是信号在时域与频域的关系，用FFT计算得到的线性卷积。
changdu1=long_N;
Blsignal3=[signal1,zeros(1,changdu1-long_N)];
Blsignal4=[signal2,zeros(1,changdu1-long_M)];
Hn1=zeros(changdu1,changdu1);
for L1=1:1:changdu1
if(L1==1)
Hn1(L1,:)=Blsignal4;
else
Hn1(L1,:)=circshift(Blsignal4,[0,1]);
Blsignal4=Hn1(L1,:);
end
end
circj=Hn1'*Blsignal3';
subplot(413)
stem(circj);
title('用函数实现的循环卷积条件为L<M+N-1');
基于MATLAB离散卷积的实现
1离散卷积实现的目的和方法
在数字信号处理(DSP)系统中，无论在时域或频域都离不开卷积运算和快速傅里叶(FFT)运算。可以使用线性卷积、循环卷积实现离散卷积，也可以根据时域信号函数与频谱函数的关系，使用快速傅里叶变换(FFT，IFFT)实现离散卷积。线性卷积是求离散系统响应的主要方法之一，许多重要应用都建立在这一理论基础上，如卷积滤波等，专用函数conv (x , h)可完成线性卷积过程。循环卷积的计算速度要远远快于线性卷积，如果选择圆周卷积的长度L≥N1+N2-1，则可以用圆周卷积取代线性卷积。FFT运算具有速度快、计算量小等优越性，是DSP的核心算法。根据信号在时域与频域的关系，同样可以用FFT计算线性卷积。通过编程实现的卷积函数深入的了解了卷积的原理。
subplot(412)
stem(signal2);
title('原信号2');
xlabel('n');ylabel('y(n)');
axis([0 18 -5 5]);
grid on;
%%%%%%%%%线性卷积函数%%%%%%%%%%%%%
JuanJi=conv(signal1,signal2);%调用函数conv实现卷积
Blsignal2=[signal2,zeros(1,changdu-long_M)];
Hn=zeros(changdu,changdu);%初始化一个矩阵
for L1=1:1:changdu
if(L1==1)
Hn(L1,:)=Blsignal2; %当处于第一行时不用进行移位操作
else
Hn(L1,:)=circshift(Blsignal2,[0,1]); %从第二行开始向右进行循环移位操作
if ((k-long_N+i>long_M)) %判断翻转序列是否超过自身的长度
fk(k)=a;
else
fk(k)=a+signal2(k-long_N+i)*signal1(long_N-i+1);
a=fk(k);
end
end
end
end
end
subplot(414)
stem(fk);
title('重新卷积函数的实现');