重叠卷积法的任意输入点数的滤波器实现.docx
《重叠卷积法的任意输入点数的滤波器实现.docx》由会员分享，可在线阅读，更多相关《重叠卷积法的任意输入点数的滤波器实现.docx（13页珍藏版）》请在冰豆网上搜索。
重叠卷积法的任意输入点数的滤波器实现
现代通信原理
设计报告
基于分段卷积法的数字滤波器
实现与C编程
实现的可行性分析
组长：
孙建东
组员：
赵玉兰
2013年12月8号
重叠卷积法的任意输入点数的滤波器实现
一、设计目的
由于信号经过低通滤波器,其本质就是实现输入信号x（n）与h（n）的积，如果x（n）序列过长则用于存储输入数据的空间很大且运算过程也相当长，且所设计的滤波器灵活性和可移植性很差,现通过研究线性卷积发现，若序列过长可采取分段卷积的思想来实现,输入长序列x（n）的滤波器算法,.，并且是一种C语言可实现、节省存储空间算法。
二、算法设计思路
首先对于四阶的数字滤波器,输入12点的序列进行卷积
x={x
（1）,x
（2）,x（3）,x（4）,x（5）,x（6,x（7）,x（8）,x（9）,x（10）,x（11）,x（12）}
h={h
（1）,h
（2）,h（3）,h（4）}
h（n）的系数关于N/2成偶对称,且h（n）是实序列
输出y（t）
h
（1）h
（2）h（3）h（4）SUMy（n）
x
（1）h
（1）*x
（1）y
（1）
x
（2）x
（1）h
（1）*x
（2）+h
（2）*x（0）y
（2）
x（3）x
（2）x
（1）h
（1）*x（3）+h
（2）*x
（1）+h（3）*x
（1）y（3）
x（4）x（3）x
（2）x
（1）h
（1）*x（4）+h
（2）*x（3）+h（3）*x
（2）+h（4）*x
（1）y（4）
x（5）x（4）x（3）x
（2）h
（1）*x（5）+h
（2）*x（4）+h（3）*x（3）+h（4）*x
（2）y（5）
x（6）x（5）x（4）x（3）h
（1）*x（6）+h
（2）*x（5）+h（3）*x（4）+h（4）*x（3）y（6）
x（7）x（6）x（5）x（4）h
（1）*x（7）+h
（2）*x（6）+h（3）*x（5）+h（4）*x（4）y（7）
x（8）x（7）x（6）x（5）h
（1）*x（8）+h
（2）*x（7）+h（3）*x（6）+h（4）*x（5）y（8）
x（9）x（8）x（7）x（6）h
（1）*x（9）+h
（2）*x（8）+h（3）*x（7）+h（4）*x（6）y（9）
x（10）x（9）x（8）x（7）h
（1）*x（10）+h
（2）*x（9）+h（3）*x（8）+h（4）*x（7）y（10）
x（11）x（10）x（9）x（8）h
（1）*x（11）+h
（2）*x（10）+h（3）*x（9）+h（4）*x（8）y（11）
x（12）x（11）x（10）x（9）h
（1）*x（12）+h
（2）*x（11）+h（3）*x（10）+h（4）*x（9）y（12）
x（12）x（11）x（10）h
（2）*x（12）+h（3）*x（11）+h（4）*x（10）y（13）
x（12）x（11）h（3）*x（12）+h（4）*x（11）y（14）
x（12）h（4）*x（12）y（15）
若将x（n）分成三段{{x
（1）,x
（2）,x（3）,x（4）},{x（5）,x（6）,x（7）,x（8）},{x（9）,x（10）,x（11）,x（12）}};分别卷积如下
h
（1）h
（2）h（3）h（4）SUMy1（n）
x
（1）h
（1）*x
（1）y1
（1）
x
（2）x
（1）h
（1）*x
（2）+h
（2）*x（0）y1
（2）
x（3）x
（2）x
（1）h
（1）*x（3）+h
（2）*x
（1）+h（3）*x
（1）y1（3）
x（4）x（3）x
（2）x
（1）h
（1）*x（4）+h
（2）*x（3）+h（3）*x
（2）+h（4）*x
（1）y1（4）
x（4）x（3）x
（2）h
（2）*x（4）+h（3）*x（3）+h（4）*x
（2）y1（5）
x（4）x（3）h（3）*x（4）+h（4）*x（3）y1（6）
x（4）h（4）*x（4）y1（7）
其他两段同理可得
y2={y2
（1）,y2
（2）,y2（3）,y2（4）,y2（5）,y2（6）,y2（7）}
y3={y3
（1）,y3
（2）,y3（3）,y3（4）,y3（5）,y3（6）,y3（7）}
发现y={y1,y2,y3}{}表示某种运算和规律
y
（1）=y1
（1）y
（2）=y1
（2）y（3）=y1（3）
y（4）=y1（4）y（5）=y1（5）+y2
（1）y（6）=y1（6）+y2
（2）
y（7）=y1（7）+y2（3）y（8）=y2（4）y（9）=y2（5）+y3
（1）
y（10）=y2（6）+y3
（2）y（11）=y2（7）+y3（3）y（12）=y3（4）
y（13）=y3（5）y（14）=y3（6）y（15）=y3（7）;
我们不难发现分段卷积的序列部分重叠相加就可得到原序列,通过大量实验发现分段卷积的序列都满足这一规律.
假设滤波器其系数为n个,输入序列m个点,且满足m=k*n,将m分成k段,分段卷积的结果为
y1={y1
（1）,y1
（2）,…y1（k）,…y1（2*n-1）}
y2={y2
（1）,y2
（2）,…y2（k）,…y2（2*n-1）}
………………
y1={yk
（1）,y
（2）,…yk（k）,…yk（2*n-1）}
最终输出的结果y（n）
y
（1）=y1
（1）……y（n）=y1（n）
y（n+1）=y1（n+1）+y2（0）…..y（2*n-1）=y1（2*n-1）+y2（n-1）y（2*n）=y2（n）
y（2*n+1）=y2（n+1）+y3（0）……y（3*n-1）=y2（2*n-1）+y3（n-1）y（3*n）=y3（n）
…….
y（（k-1）*n+1）=yk-1（n+1）+yk（0）……y（k*n-1）=yk-1（2*n-1）+yk（n-1）
y（k*n）=yk（n）;
y（k*n=1）=yk（n+1）…….y（k*n+n-1）=yk（2*n-1）
同时又因为序列的对称性在存储滤波器系数是可知存取前n/2个点，后n/2个点通过倒序的算法得出,同样运用分段卷积的思想以h（n）={1,2,2,1}为例可以得出如下结果
1221SUMz（n）
x
（1）1*x
（1）z
（1）
x
（2）x
（1）1*x
（2）+2*x（0）z
（2）
x（3）x
（2）x
（1）1*x（3）+2*x
（1）+2*x
（1）z（3）
x（4）x（3）x
（2）x
（1）1*x（4）+2*x（3）+2*x
（2）+1*x
（1）z（4）
x（4）x（3）x
（2）2*x（4）+2*x（3）+1*x
（2）z（5）
x（4）x（3）2*x（4）+1*x（3）z（6）
x（4）1*x（4）z（7）
将h（n）={1221},分为h1（n）={12};h2（n）={21}分别卷积x（n）得出如下结果
12SUMz1（n）z（n）
x
（1）1*x
（1）z1
（1）z
（1）
x
（2）x
（1）1*x
（2）+2*x
（1）z1
（2）z
（2）
x（3）x
（2）1*x（3）+2*x
（2）z1（3）+z（3）
x（4）x（3）1*x（4）+2*x（3）z1（4）+z（4）
x（4）2*x（4）z1（5）+z（5）
21SUMz2（n）
x
（1）2*x
（1）z2
（1）
x
（2）x
（1）2*x
（2）+1*x
（1）z2
（2）
x（3）x
（2）2*x（3）+1*x
（2）z2（3）
x（4）x（3）2*x（4）+1*x（3）z2（4）z（6）
x（4）1*x（4）z2（5）z（7）
同理可得出
当输入滤波器的序列分为h1（n）和h2（n）都为n点,输入序列为m点
z（n）=h（n）*x（n）
z1（n）=h1（n）*x（n）;z2（n）=h2（n）*x（n）
z1={z1
（1）,z1
（2）,…z1（k）…,z1（n+m-1）};
z2={z2
（1）,z2
（2）,…z2（k）…,z2（n+m-1）};
z={z
（1）,z
（2）,…z（k）…,z（n+m-1）…z（2n+m-1）}
z与z1,z2的关系
z
（1）=z1
（1）,….z（n）=z1（n）
z（n+1）=z1（n+1）+z2（0）……z（n+m-1）=z1（n+m-1）+z（m-1）
z（n+m）=z2（m）……z（2*n+m-1）=z（n+m-1）
x（n）*h（n）的线性卷积的实现办法
定义离散序列线性卷积来看，就三个步骤移位、乘积、求和。
只要能很好的完成这三个步骤线性卷积就可以实现。
其实现办法，圆周卷积实现线性卷积的思想，完成n点的x（n）与m点的h（n）线性卷积，其本质就是将x（n）扩成长度为n+m-1的x1（n），h（n）扩成长度为n+m-1的h1（n）的两序列的圆周卷积得到长度为n+m-1的序列y（n）；
例y（n）=x（n）*h（n）；x（n）={x（0）,x
（1）,x
（2）,x（3）}；h（n）={h（0）,h
（1）,h
（2）,h（3）};
h1（n）000h（0）h
（1）h
（2）h（3）y（n）=h（n）*x（n）
x1（n）x（3）x
（2）x
（1）x（0）000
x1（n-1）0x（3）x
（2）x
（1）x（0）00
x1（n-2）00x（3）x
（2）x
（1）x（0）0
x1（n-3）000x（3）x
（2）x
（1）x（0）y（i）=sum（h（i）*x（i））
x1（n-4）x（0）000x（3）x
（2）x
（1）i=0,1,2,3,4,5,6
x1（n-5）x
（1）x（0）000x（3）x
（2）
x1（n-6）x
（2）x
（1）x（0）000x（3）
三、该算法的编程实现
基于以上的算法设计出如下程序
程序一：
完成圆周移位函数名：
clc_shift
%函数名称:
clc_shift完成对x的循环移位
%x:
需要移位的序列n:
移位序列的长度
functionout=clc_shift（x,n）
temp=x（1,n）;%将数组最后一位的元素存入暂存器
fori=1:
n-1
x（1,n-i+1）=x（1,n-i）;%数组中所有元素向后移一位
end
x（1,1）=temp;%将数组最后一个元素赋给第一位
out=x;%将移位后的序列输出
end
程序二：
完成线性卷积函数名：
conv_user
%函数名称:
conv_user完成x与h的线性卷积
%x,h:
参与线性卷积的序列n:
序列x的长度m:
序列h的长度
functiony=conv_user（x,h,n,m）
x1=zeros（1,n+m-1）;%开辟长度为n+m-1的序列x1
h1=zeros（1,n+m-1）;%开辟长度为n+m-1的序列h1
y=zeros（1,n+m-1）;%开辟长度为n+m-1的序列y
fori=1:
n+m-1%给x1初始化数据
ifi<=n
x1（1,i）=x（1,n-i+1）;%前n个点将x的序列反向存储
else
x1（1,i）=0;%后m-1个点赋0
end
end
%x1=[x2,zeros（1,m-1）];
forj=1:
m%给h1初始化数据
h1（1,n+j-1）=h（1,j）;%前n-1点赋0，后m点顺序赋h的值
end
fori=1:
n+m-1%x1序列移位n+m-1次
sum=0;
forj=1:
n+m-1%完成一次移位得到y（i）的值
sum=sum+x1（1,j）*h1（1,j）;%y（i）的各项乘积求和
end
y（1,i）=sum;
x1=clc_shift（x1,n+m-1）;%完成一次运算该数组移位一次
end
end
程序三:
实现对称卷积,函数名:
y=juanji（h,x）
Functiony=juanji（h,x）%实现偶对称序列的一种卷积算法
m=length（x）;%输入x序列的长度
n=length（h）;%系统函数的序列前n/2的长度
y=zeros（1,2*n+m-1）;%开辟输出序列的长度
y1=[0:
n+m-1];y2=[0:
n+m-1];%开辟中间序列的数据长度
y1=conv_user（x,h）;%实现h（n）的前n/2点的卷积
fori=0:
n/2-1%转换成后关于n/2偶对称后h（n）后n/2的序列
t=h（i+1）;
h（i+1）=h（n-i）;%实现前后数据交换
h（n-i）=t;
end
y2=conv_user（x,h）;%实现后n/2个点的卷积
%以下实现偶对称的两个子序列卷积数据通过某种规律
forii=1:
n%前输出序列n个点
y（ii）=y1（ii）;
end
forii=n+1:
n+m-1%输出序列n+1个点到n+m-2个点
y（ii）=y1（ii）+y2（ii-n）;%原序列与两个子序列的关系
end
forii=n+m:
2*n+m-1%输出序列最后n个点
y（ii）=y2（ii-n）;
end
end
程序四:
实现输入序列x（n）与h（n）的卷积,函数名:
functiony=FIR_L（x,h）
%设计n阶FIR的低通滤波器输入点数N为任意点数
functiony=FIR_L（x,h）%x输入采样点的向量,y输出序列,h为所设计低通滤波器的前0.5*n个滤波器系数
m=length（x）;%统计x的序列长度
n=2*length（h）;%确定h的长度即滤波器的阶数
y=zeros（1,n+m-1）;%y输出序列的长度
k=ceil（m/n）;%确定m与n的倍数k
xi=zeros（1,k*n）;
xi=[x,zeros（1,k*n-m）];
x0=zeros（k,n）;%x输入序列分k段存储在x0
y0=zeros（k,2*n-1）;%y0存储分段卷积的结果
forjj=1:
k
x0（jj,:
）=xi（1,（jj-1）*n+1:
jj*n）;%将x序列分成k段
y0（jj,:
）=juanji（h,x0（jj,:
））;%实现分段卷积
end
forjj=1:
k
ifjj==1
forii=1:
n
y（ii）=y0（jj,ii）;%第一段的存储
end
elseif1forii=1:
n-1
y（（jj-1）*n+ii）=y0（（jj-1）,（n+ii））+y0（jj,ii）;%中间段的存储
end
y（jj*n）=y0（jj,n）;
end
ifjj==k
forii=1:
n-1
y（jj*n+ii）=y0（jj,ii+n）;%最后一段的存储
end
end
end
end
四、仿真验证
利用MATLAB工具设计一个输入点数任意的20阶数字低通滤波器要求
Wp=40Hz,Ws=60Hz,Fs=200Hz,现有一信号x=x0+x1;x0=cos（2*pi*20*t）是有用信号x1是随机干扰信号,通过该滤波器后滤除干扰信号保留有用信号,保证通带和阻带内纹波较小
滤波器h（n）设计
利用FDATOOL获得滤波器的h（n）如图
编写验证程序
输入序列为100点的序列采样频率Fs=200Hz,采样时间0.5s调用MATKLAB内部函数编写仿真M文件
t0=[0:
1/200:
99/200];%输入序列的时间长度
t1=[0:
1/200:
118/200];%输出序列的时间长度
x0=cos（2*pi*20*t0）;%有用输出信号
x1=sin（2*pi*60*t0）;%无用干扰信号
x=x0+x1;%输入信号
h1=[0.00610.0121-0.0155-0.01910.02800.0385-0.0557-0.08290.14650.44830.44830.1465-0.0829-0.05570.03850.0280-0.0191-0.01550.01210.0061];%20阶低通滤波器系数
y=conv（x,h1）;%信号过低通滤波器即实现输入信号与输出信号的卷积
figure;%
subplot（2,2,1）;
plot（t0,x0,'r','LineWidth',2）;%显示有用信号
axis（[0.1,0.4,-2,2]）;
title（'x0的波形'）;holdon
subplot（2,2,2）;
plot（t0,x1,'b','LineWidth',2）;%显示无用信号
axis（[0.1,0.3,-2,2]）;
title（'x1的波形'）;holdon
subplot（2,2,3）;
plot（t0,x,'g','LineWidth',2）;%显示输入信号
axis（[0.1,0.4,-2,2]）;
title（'x的波形'）;holdon
subplot（2,2,4）;
plot（t1,y,'r','LineWidth',2）;%显示输出信号
axis（[0.1,0.4,-2,2]）;
title（'y的波形'）;holdon
程序运行结果:
程序五:
测试该滤波器的性能,函数名称:
FIRtest
t0=[0:
1/200:
99/200];%输入序列仿真时间长度
x0=cos（2*pi*20*t0）;%有用输出信号
subplot（2,2,1）;
plot（t0,x0,'r'）;%显示有用信号
axis（[0.1,0.4,-2,2]）;title（'x0的波形'）;
x1=sin（2*pi*60*t0）;%无用干扰信号
subplot（2,2,2）;
plot（t0,x1,'y'）;%显示无用信号
axis（[0.1,0.3,-2,2]）;title（'x1的波形'）;
x=x0+x1;%输入信号的函数表达式
subplot（2,2,3）;
plot（t0,x,'b'）;%显示输入序列
axis（[0.1,0.5,-2,2]）;title（'输入序列x（n）'）;
t0=[0:
1/200:
118/200];
h=[-0.05090.00860.04150.04830.0051-0.0550-0.05780.04440.21460.3461];%所设计滤波器的h（n）的系数向量
y=FIR_L（x,h）;%输入信号通过滤波器实现滤波,即x与y的卷积
subplot（2,2,4）;
plot（t0,y,'r'）;%显示滤波后的输出序列
axis（[0.1,0.4,-2,2]）;title（'滤波输出序列y（n）'）;
仿真结果如图:
六、仿真结果分析
通过对该算法的仿真有效的验证，该算法用C语言实现具有一定可行性，该算法的优势在于该滤波器对输入信号为长序列实现卷积的结构固定但输入序列长度可变只要输入序列满足输入点数m=kn，n为滤波器的阶数，都可以实现有效的卷积；且分段卷积理论上可以实现与无限长输入序列卷积，但实际受到实际器件资源的限制无限长是做不到的，在C语言实现过程中所有接口采用宏定义当输入序列变化时只需改变相应宏定义的参数值，不必改动程序的结构，这样在程序的可移植性方面有很好的优势。
上述所编写的算法在存储空间的分配上还可以进一步的优化，由于时间原因本人这次仅为了验证该算法的可行性，而并未过多思考如何进一步优化该算法对输入序列存储空间分配的优化，今后可进一步研究改进该算法。
升级会员 