CN104331862A - 一种并联型分数阶零相位滤波器及其滤波方法 - Google Patents

Abstract

本发明公开了一种并联型分数阶零相位滤波器及其滤波方法，属于数字信号处理技术领域。包括前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器以及比例运算器；所述输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器进行滤波，滤波后的信号通过比例运算器调节后得到输出信号。本发明计算简便，应用于数字信号处理领域，滤波结果无相位失真，在信号去噪和信息保留之间达到有效的平衡，并且可以依据实际需要选择阶次。

Description

一种并联型分数阶零相位滤波器及其滤波方法

技术领域

本发明涉及一种并联型分数阶零相位滤波器及其滤波方法，属于数字信号处理技术领域。

背景技术

信号去噪处理是信号处理中的一个重要并且基础的环节。到目前为止，学者们提出了诸多去噪算法，包括均值滤波、顺序统计滤波、低通滤波和维纳滤波等等滤波算法。上述算法可以在不同程度上降低噪声，但是在滤波过程中却存在丢失信号中部分有用信息的风险，从而导致信号模糊。产生这一现象的原因是上述去噪算法都直接或间接的将其去噪模型建立在整数阶次积分上，而整数阶次的积分对于高频信息的抑制非常明显，当信号的有效频带与噪声的频带产生重叠的时候，常会导致信号误滤波。因此，提出一种同时具备抑制噪声和保留信号有用成分的滤波器方法迫在眉睫。

近三百年来，分数阶微积分在数学分析领域中已成为一个重要分支，但对于大多数工程技术界学者而言还鲜为人知。图像处理领域中，学者们近年来相继提出诸多基于分数阶微积分的图像处理算法，但是在一维信号处理领域，尤其是生物信号处理和非因果系统领域，如批信号处理，分数阶微积分的应用仍是一个急需研究的新兴学科分支。

分数阶积分运算在对一维信号进行处理的时候既能非线性加强信号中的低频和中频成分，抑制高频成分，又能在一定程度上尽量多地保留信号的高频成分。分数阶运算的定义并不唯一。欧氏空间下最常使用的定义是Riemann-Liouville和Grünwald-Letnikov定义。由于Grünwald-Letnikov定义基于有限数量的离散的点，而数字信号是由离散值组成的，Grünwald-Letnikov定义被广泛应用在数字信号处理中。

发明内容

发明目的：为了克服现有技术中存在的不足，本发明提供一种并联型分数阶零相位滤波器及其滤波方法，该滤波器及其滤波方法能够解决一维信号在滤波过程中滤除噪声的同时将与噪声频段重叠的有用信号滤除的问题，其滤波结果无相位失真。

为实现上述目的，本发明采用的技术方案为：一种并联型分数阶零相位滤波器，包括中心Grünwald-Letnikov微积分算子滤波器所述中心Grünwald-Letnikov微积分算子滤波器包括前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器以及比例运算器，所述前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器相互并联在一起之后再与比例运算器串联；所述输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器进行滤波，滤波后的信号通过比例运算器调节后得到输出信号。

优选的：所述中心Grünwald-Letnikov微积分算子滤波器与前向Grünwald-Letnikov微积分算子滤波器和后向Grünwald-Letnikov微积分算子滤波器的关系为：

$D_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\π}}{2}\right)}(D_l^v+D_r^v),v\≠2n+1,n\∈Z$

其中v为微积分阶次，取任意非奇实数。

优选的：所述输出信号为：

Y(e^iω)＝X(e^iω)·|ω|^v

其中v为分数阶微积分滤波器的阶次，取任意非奇实数；X(e^iω)为输入信号；i为虚数单位，ω为信号频率。

优选的：所述中心Grünwald-Letnikov微积分算子滤波器前向Grünwald-Letnikov微积分算子滤波器以及后向Grünwald-Letnikov微积分算子滤波器均通过差分方程的形式来构建。

一种并联型分数阶零相位滤波器的滤波方法，首先输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器进行滤波；然后将滤波后的信号通过比例运算器调节后得到输出信号。

优选的：所述中心Grünwald-Letnikov微积分算子滤波器与前向Grünwald-Letnikov微积分算子滤波器和后向Grünwald-Letnikov微积分算子滤波器的关系为：

$D_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\π}}{2}\right)}(D_l^v+D_r^v),v\≠2n+1,n\∈Z$

其中v为微积分阶次，取任意非奇实数。

优选的：所述输出信号为：

Y(e^iω)＝X(e^iω)·|ω|^v

其中v为分数阶微积分滤波器的阶次，取任意非奇实数；X(e^iω)为输入信号；i为虚数单位，ω为信号频率。

优选的：所述中心Grünwald-Letnikov微积分算子滤波器前向Grünwald-Letnikov微积分算子滤波器以及后向Grünwald-Letnikov微积分算子滤波器均通过差分方程的形式来构建。

本发明提供的并联型分数阶零相位滤波器及其滤波方法，相比现有技术，具有以下有益效果：

(1)在一维信号处理过程中滤除噪声的同时，可以更好地保留与噪声频段重叠的信号有用成分。在信号去噪和信息保留之间达到有效的平衡。

(2)设计方法简便，通过对信号进行分数阶微积分运算，分别实现对信号的分数阶前向滤波和后向滤波，通过并联设计最终得到的分数阶零相位滤波器设计方法简便且算法效率高。

(3)由于信号通过前向Grünwald-Letnikov微积分算子相角超前而信号通过后向Grünwald-Letnikov微积分算子相角滞后因此本发明无失真的还原信号特征的位置。

(4)设计灵活性大，通过调节分数阶微积分的阶次，可以依据实际需要的设计指标选择合适的阶次。

综上所述：将本发明基于前向Grünwald-Letnikov微积分算子和后向Grünwald-Letnikov微积分算子的并联型分数阶零相位滤波器应用于数字信号处理领域，可达到零相位滤波效果。同时本发明计算简便，涉及分数阶微积分运算；滤波结果无相位失真；在信号去噪和信息保留之间达到有效的平衡，并且可以依据实际需要选择阶次。本发明的分数阶零相位滤波器可应用于脑电、肌电、眼电等生物电信号以及其他非平稳信号处理，尤其适用于批信号处理领域。

附图说明

图1为本发明并联型分数阶零相位滤波器的流程图。

图2为本发明并联型分数阶零相位滤波器的信号处理流程图。

图3为卷积模板曲线图，其中v分别取-0.5和0.5两种情况，m＝100。

具体实施方式

下面结合附图对本发明作更进一步的说明。

一种并联型分数阶零相位滤波器，如图1、2所示，包括中心Grünwald-Letnikov微积分算子滤波器所述中心Grünwald-Letnikov微积分算子滤波器包括前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器以及比例运算器，所述前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器相互并联在一起之后再与比例运算器串联；所述输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器和后向Grünwald-Letnikov微积分算子滤波器进行滤波，滤波后的信号通过比例运算器调节后得到输出信号。

如图1、2所示，本发明设计的分数阶零相位滤波器，将前向Grünwald-Letnikov微积分算子和后向Grünwald-Letnikov微积分算子并联连接，再通过比例运算，得到基于中心Grünwald-Letnikov微积分算子的零相位滤波器。输入信号通过基于中心Grünwald-Letnikov微积分算子的零相位滤波器后，中间输出信号再次通过比例环节，得到零相位失真的输出信号。

本发明所解决的技术问题在于获得一种比现有的基于整数阶积分的滤波方法更理想的信号处理算法，其可以在信号处理过程中杜绝相位失真，并在滤除噪声的同时，保留信号中与噪声频段重叠的有用成分。在具体说明本发明内容之前，有必要对本说明书所用符号涵义进行简要说明。

根据一阶左导数的定义，可得到:

$D_{l}f\left(t\right)=\underset{h\→0}{\lim }\frac{f\left(t\right)-f(t-h)}{h}$

其中，f(t)为输入信号。

将定义为一阶n重左导数，可得到:

$D_l^{n}f\left(t\right)=\underset{h\→0}{\lim }\frac{1}{h^n}\underset{k=0}{\overset{n}{\mathrm{\Σ}}}{(-1)}^k\left(\begin{array}{c}n\\ k\end{array}\right)f(t-\mathrm{kh})$

其中， $\left(\begin{array}{c}n\\ k\end{array}\right)=\frac{n!}{k!(n-k)!}$ 为二项式系数。

当n取分数时，该二项式系数可替换为 $\left(\begin{array}{c}v\\ k\end{array}\right)=\frac{v(v-1)(v-2)\·\·\·(v-k+1)}{k!},$ 其中v为任意非整数。由此我们定义v阶前向Grünwald-Letnikov微积分算子用符号表示。对函数f(t)进行v阶前向Grünwald-Letnikov微积分，得到下式：

$D_l^{v}f\left(t\right)=\underset{h\→0}{\lim }\frac{1}{h^n}\underset{k=0}{\overset{\left[\frac{t-a}{h}\right]}{\mathrm{\Σ}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)f(t-\mathrm{kh})$

其中，t和a分别是前向Grünwald-Letnikov微积分算子的上限和下限。该微积分算子将n阶微分与n阶积分运算统一用一个微积分算子来表示，且将整数域的运算推广到非整数域。当v>0时，上式表示微分运算，当v<0时，上式表示积分运算。微积分算子由此得名。信号进行v阶前向Grünwald-Letnikov微积分的傅里叶变换表示为:

$F\left[D_l^{v}f\left(t\right)\right]={\left(\mathrm{i\&omega;}\right)}^{v}f\left(\mathrm{\&omega;}\right)$

因此，v阶前向Grünwald-Letnikov微积分算子的频域函数为:

$H_l^v\left(e^{\mathrm{iw}}\right)={\left(\mathrm{i\&omega;}\right)}^v={\left|\mathrm{\&omega;}\right|}^v{e}^{i\frac{\mathrm{v\&pi;}}{2}\mathrm{sgn}\left(\mathrm{\&omega;}\right)}$

对应的幅频和相频响应分别为：

根据一阶右导数的定义，可得到

$D_{r}f\left(t\right)=\underset{h\&RightArrow;0}{\lim }\frac{f\left(t\right)-f(t+h)}{h}$

其中f(t)为输入信号。

将定义为一阶n重右导数，可得到:

$D_r^{v}f\left(t\right)=\underset{h\&RightArrow;0}{\lim }\frac{1}{h^v}\underset{k=0}{\overset{\left[\frac{b-t}{h}\right]}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)f(t-\mathrm{kh})$

其中, $\left(\begin{array}{c}n\\ k\end{array}\right)=\frac{n!}{k!(n-k)!}$ 为二项式系数。

当n取分数时，该二项式系数可替换为 $\left(\begin{array}{c}v\\ k\end{array}\right)=\frac{v(v-1)(v-2)\&CenterDot;\&CenterDot;\&CenterDot;(v-k+1)}{k!},$ 其中v为任意非整数。由此我们定义v阶后向Grünwald-Letnikov微积分算子用符号表示。对函数f(t)进行v阶后向Grünwald-Letnikov微积分，得到下式：

$D_r^{v}f\left(t\right)=\underset{h\&RightArrow;0}{\lim }\frac{1}{h^v}\underset{k=0}{\overset{\left[\frac{b-t}{h}\right]}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)f(t-\mathrm{kh})$

其中，b和t分别是后向Grünwald-Letnikov微积分算子的上限和下限。信号进行v阶后向Grünwald-Letnikov微积分的傅里叶变换表示为:

$F\left[D_r^{v}f\left(t\right)\right]={(-\mathrm{i\&omega;})}^{v}f\left(\mathrm{\&omega;}\right)$

因此，v阶后向Grünwald-Letnikov微积分算子的频域函数为:

$H_r^v\left(e^{\mathrm{iw}}\right)={(-\mathrm{i\&omega;})}^v={\left|\mathrm{\&omega;}\right|}^v{e}^{-i\frac{\mathrm{v\&pi;}}{2}\mathrm{sgn}\left(\mathrm{\&omega;}\right)}$

对应的幅频和相频响应分别为：

前向Grünwald-Letnikov微积分算子和后向Grünwald-Letnikov微积分算子可被视为两种特殊的滤波器，其幅值响应随信号频率和分数阶微积分阶次变化而变化，信号通过前向Grünwald-Letnikov微积分算子相角超前而信号通过后向Grünwald-Letnikov微积分算子相角滞后

本发明设计的分数阶零相位滤波器，将前向Grünwald-Letnikov微积分算子和后向Grünwald-Letnikov微积分算子并联连接，通过比例运算，得到基于中心Grünwald-Letnikov微积分算子的零相位滤波器。输入信号通过基于中心Grünwald-Letnikov微积分算子的零相位滤波器后，中间输出信号再次通过比例环节得到零相位失真的输出信号。用数学表达式来说明与的关系，如下式：

$D_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}(D_l^v+D_r^v),v\&NotEqual;2n+1,n\&Element;Z$

因此，对函数f(t)进行v阶中心Grünwald-Letnikov微积分可以表示为：

$D_c^{v}f\left(t\right)=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}\underset{h\&RightArrow;0}{\lim }\frac{1}{h^v}\underset{k=0}{\overset{N}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)[f(t-\mathrm{kh})+f(t+\mathrm{kh})]$

信号f(t)进行v阶中心Grünwald-Letnikov微积分的傅里叶变换为:

$F\left[D_c^{v}f\left(t\right)\right]=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}[{\left(\mathrm{i\&omega;}\right)}^v+{(-\mathrm{i\&omega;})}^v]={\left|\mathrm{\&omega;}\right|}^v$

因此，v阶中心Grünwald-Letnikov微积分算子的频域函数为:

$H_c^v\left(e^{\mathrm{i\&omega;}}\right)=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}[{\left(\mathrm{i\&omega;}\right)}^v+{(-\mathrm{i\&omega;})}^v]={\left|\mathrm{\&omega;}\right|}^v$

对应的幅频和相频响应分别为：

具体实施时，由于数字信号处理是基于离散信号进行的，因此我们需将上述时域函数进行离散化操作，即用差分方程的形式来构建数字信号滤波器。取h＝1，且当或足够大的时候，可省略公式中的取极限这一步骤。前向Grünwald-Letnikov微积分算子后向Grünwald-Letnikov微积分算子和中心Grünwald-Letnikov微积分算子的差分近似定义如下：

${\mathrm{\&Delta;}}_l^{v}f\left(t\right)=\underset{k=0}{\overset{N}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)f(t-k)$

${\mathrm{\&Delta;}}_r^{v}f\left(t\right)=\underset{k=0}{\overset{N}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)f(t+k)$

${\mathrm{\&Delta;}}_c^{v}f\left(t\right)=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}\underset{k=0}{\overset{N}{\mathrm{\&Sigma;}}}{(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right)[f(t-k)+f(t+k)],$ 其中v≠2n+1,n∈Ζ

令 ${(-1)}^k\left(\begin{array}{c}v\\ k\end{array}\right){(-1)}^k\frac{v(v-1)\&CenterDot;\&CenterDot;\&CenterDot;(v-k+1)}{k!}\stackrel{\mathrm{\&Delta;}}{=}{a}_k$

为了计算出前向Grünwald-Letnikov微积分算子的差分方程，取卷积模板 $M_l^v=[a_m\&CenterDot;\&CenterDot;\&CenterDot;a_2{a}_{1}1],$ 差分方程为 ${\mathrm{\&Delta;}}_l^{v}f\left(t\right)=f\left(t\right)*M_l^v.$

同理，为了计算出后向Grünwald-Letnikov微积分算子的差分方程，取卷积模板 $M_r^v=[1a_1{a}_2\&CenterDot;\&CenterDot;\&CenterDot;a_m],$ 差分方程为 ${\mathrm{\&Delta;}}_r^{v}f\left(t\right)=f\left(t\right)*M_r^v.$

当中心Grünwald-Letnikov微积分算子应用在有限长度的信号上时

${\mathrm{\&Delta;}}_c^{v}f\left(t\right)=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}[{\mathrm{\&Delta;}}_l^{v}f\left(t\right)+{\mathrm{\&Delta;}}_r^{v}f\left(t\right)]$

或者

${\mathrm{\&Delta;}}_c^{v}f\left(t\right)f\left(t\right)*\left\{\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}\right[M_l^v\left(t\right)+M_r^v\left(t\right)\left]\right\}$

或者

${\mathrm{\&Delta;}}_c^{v}f\left(t\right)=f\left(t\right)*M_c^v\left(t\right)$

其中， $M_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}[a_m\&CenterDot;\&CenterDot;\&CenterDot;a_2{a}_{1}2a_1{a}_2\&CenterDot;\&CenterDot;a_m]$

图3所示为卷积模板曲线图，其中v分别取-0.5和0.5两种情况，m＝100。图3左半边表示右半边表示两条卷积模板曲线偶对称且在两端趋于0，表现了零相位滤波器的通性，且保证了适度的计算复杂度。

所述分数阶零相位滤波器的输出信号采用如下公式表示：

Y(e^iω)＝X(e^iω)·|ω|^v，其中，v为分数阶微积分滤波器的阶次，取任意非奇实数；X(e^iω)为输入信号；i为虚数单位，ω为信号频率。

本发明涉及的微积分的阶次不是传统的整数阶，而是非整数阶，工程应用中一般取分数或者有理小数。本发明提出的分数阶零相位滤波器设计方法结构简单，运算规则简明，计算量小，效率高。将该发明技术应用于心电信号工频干扰去噪、肌电干扰去噪以及合成干扰去噪均取得了较好的效果，在抑制高频噪声以及波形保持方面优于其它方法。此外，本项目所设计的分数阶零相位滤波器还可广泛应用于脑电、肌电、眼电等生物电信号以及其他非平稳信号处理，该研究成果必将提升我国在该领域的科研学术和实际应用水平。

以上所述仅是本发明的优选实施方式，应当指出：对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (8)

1.一种并联型分数阶零相位滤波器，其特征在于：包括中心Grünwald-Letnikov微积分算子滤波器所述中心Grünwald-Letnikov微积分算子滤波器包括前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器以及比例运算器，所述前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器相互并联在一起之后再与比例运算器串联；所述输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器进行滤波，滤波后的信号通过比例运算器调节后得到输出信号。

2.根据权利要求1所述的并联型分数阶零相位滤波器，其特征在于：所述中心Grünwald-Letnikov微积分算子滤波器与前向Grünwald-Letnikov微积分算子滤波器和后向Grünwald-Letnikov微积分算子滤波器的关系为：

$D_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}(D_l^v+D_r^v),v\&NotEqual;2n+1,n\&Element;Z$

其中，v为微积分阶次，取任意非奇实数。

3.根据权利要求2所述的并联型分数阶零相位滤波器，其特征在于：所述输出信号为：

Y(e^iω)＝X(e^iω)·|ω|^v

其中，v为分数阶微积分滤波器的阶次，取任意非奇实数；X(e^iω)为输入信号；i为虚数单位，ω为信号频率。

4.根据权利要求3所述的并联型分数阶零相位滤波器，其特征在于：所述中心Grünwald-Letnikov微积分算子滤波器前向Grünwald-Letnikov微积分算子滤波器以及后向Grünwald-Letnikov微积分算子滤波器均通过差分方程的形式来构建。

5.一种基于权利要求1所述的并联型分数阶零相位滤波器的滤波方法，其特征在于：首先输入信号分别通过前向Grünwald-Letnikov微积分算子滤波器后向Grünwald-Letnikov微积分算子滤波器进行滤波；然后将滤波后的信号通过比例运算器调节后得到输出信号。

6.根据权利要求5所述的并联型分数阶零相位滤波器的滤波方法，其特征在于：所述中心Grünwald-Letnikov微积分算子滤波器与前向Grünwald-Letnikov微积分算子滤波器和后向Grünwald-Letnikov微积分算子滤波器的关系为：

$D_c^v=\frac{1}{2\cos \left(\frac{\mathrm{v\&pi;}}{2}\right)}(D_l^v+D_r^v),v\&NotEqual;2n+1,n\&Element;Z$

其中v为微积分阶次，取任意非奇实数。

7.根据权利要求6所述的并联型分数阶零相位滤波器的滤波方法，其特征在于：所述输出信号为：

Y(e^iω)＝X(e^iω)·|ω|^v

其中v为分数阶微积分滤波器的阶次，取任意非奇实数；X(e^iω)为输入信号；i为虚数单位，ω为信号频率。

8.根据权利要求7所述的并联型分数阶零相位滤波器的滤波方法，其特征在于：所述中心Grünwald-Letnikov微积分算子滤波器前向Grünwald-Letnikov微积分算子滤波器以及后向Grünwald-Letnikov微积分算子滤波器均通过差分方程的形式来构建。