沈政伟等：分数阶对称平移不变过完备小波的构造及其在轴承故障诊断中的应用 381· 则A(x)可以表示为 构条件的等价多相位矩阵表达式： A(x）=- 1-B(x)x* (16) [H(z) Gg81-4 (20) 式中 A(x)可以通过泰勒级数展开式得到： rHo(z） Ho(2)Ho.(2)1 =A(x)+0(1x).(17) H(）= 1-)(-) [Ho(2) H,(e)H2(a)] TGm(z） Go(=)Go (= 展开式子 G(z)= Gio(z) Gu(2) G2(z） LGm(a)G(z)G2(a)」 通过式(20)得到 和 G"(a-)G(z)=I,-H(2)H(a). (21) 下面将通过矩阵谱因式分解的方式求解上式中的 G(z).对于矩阵谱因式分解，文献3]给出了多种方 此处的 法的详细介绍，但这些方法仅限于简单例子的手工计 N-1 ]为组合运算，取泰勒展开式中的前 算.文献2]中使用的是对称因式分解方法，但计算结 K项，即可得到A(x).再通过标准化(H(1)=√6求 果还需要最后通过仿酉矩阵调整才能使高通滤波器组 出T(x),最后经过参数转换即得到如下结果： 满足近似平移不变特性，因此计算复杂度非常高.本 文将采用拓普利兹矩阵分解法对式(21)进行直接分 解，该方法是一种迭代计算方法，计算方法非常简单. q*2 -g+N-1 通过程序计算对比发现，本文采用拓普利兹矩阵 -1 N-1 分解法求高通滤波器组具有以下三点优势：(1)对称 因式分解法回和使用拓普利兹矩阵分解法迭代5000 (告)”(2 (18) 次计算时间开销虽然接近，但是本文方法在经过3000 次拓普利兹矩阵分解迭代运算（时间为5000次迭代 式中，z-”因子是用来调整平移量.同理当V为偶 的1/2)之后，其计算结果误差范围已经在0.0055%~ 数时，也可以得到相应的表达式： 0.0097%.(2)利用拓普利兹矩阵分解法得到的高通 滤波器组已经满足近似平移不变特性，无需再如文献 ]一样做额外的计算，因此省去繁杂的调整过程，降 +22 低了计算复杂度.(3)根据不同的数据精度需求和滤 p-q+N-1Y N-1 波器的设计复杂程度，采用拓普利兹矩阵分解法也具 2-2 有更好的灵活度，同时也更方便计算机实现. 表1是文献2]N=3和K=1高通滤波器组.对 (19) 比表2中由本文中提出的分解法经过3000次迭代运 因此给定N和K,通过式(18)和式(19)可以很容 算之后的结果，可以看出两种方法求出的高通滤波器 易求出相应的低通滤波器。需要强调说明的是，利用 组是近似一样的.这验证了本文方法的正确性. 式(18)和式(19)构造的分数阶低通滤波器是具有线 表1文献2]中的N=3和K=1的高通滤波器系数 性相位特性的（对称性）. Table 1 High-pass filter coefficients of N=3 and K=1 in Ref. 2.2高通滤波器的构造 50(n) &(n) &2(n) 构造了具有线性相位特性的低通滤波器之后，将 0.64917778505741 0 求解相应的高通滤波器，从而完成图2中滤波器组的 0 -0.482626543662260.63770868747435 0 构造.首先，将低通滤波器H(z)用多相位方法表示为 H(z)=H。()+z3H,(),则低通滤波器部分被分解 -0.15059130119969-0.466878032128120.64520631583316 成两个低通通道H。(z)和H,(z).为方便使用完全可 -0.01477135217528-0.14815175304968-0.49098922627425 重构条件(2)~(4)，再将H。(z)展开成多相位形式， -0.00118858802016-0.02267890229656-0.13149490989732 即H(z)=Ho(z3)+zHo(z3)+z2H2(z),H,(z）和 0 0 -0.02152627546889 G,(z)(i=0,1,2)按照同样的方式分解，得到完全可重 0 0 -0.00119590419272沈政伟等: 分数阶对称平移不变过完备小波的构造及其在轴承故障诊断中的应用 则 A( x) 可以表示为 A( x) = 1 － B( x) xK ( 1 － x) N/ ( 2 3 4 － ) x N ． ( 16) A( x) 可以通过泰勒级数展开式得到: 1 ( 1 － x) N/ ( 2 3 4 － ) x N = A( x) + O( | x | K ) ． ( 17) 展开式子 1 ( 1 － x) N = ∑ ∞ n = 0 n + N － 1 ( N ) － 1 xn 和 1 ( 3 /4 － x) N = ∑ ∞ n = 0 ( 3 /4) N n + N － 1 ( N ) － 1 ( 3 /4) n xn ， 此处的 n + N － 1 ( N ) － 1 为组合运算，取泰勒展开式中的前 K 项，即可得到 A( x) ． 再通过标准化( H( 1) = 6槡) 求 出 T( x) ，最后经过参数转换即得到如下结果: H( z) = z － ( K － 1) 槡 ( 6 1 + z － 1 ) 2 ( N 1 + z － 1 + z － 2 ) 3 N · ∑ K －1 p = 0 ∑ p q = 0 q + N 2 － 1 N 2           － 1 p － q + N － 1 ( N ) － 1 ( · ) 4 3 p － ( q － z － 1 + 2 － z ) 4 p ． ( 18) 式中，z － ( K － 1) 因子是用来调整平移量． 同理当 N 为偶 数时，也可以得到相应的表达式: H( z) = z － ( K － 1) 槡 ( 6 1 + z － 1 ) 2 N ( + 1 1 + z － 1 + z － 2 ) 3 N · ∑ K －1 p = 0 ∑ p q = 0 q + N 2 － 1 2 N 2 －           1 2 p － q + N － 1 ( N ) － 1 ( · ) 4 3 p － ( q － z － 1 + 2 － z ) 4 p ． ( 19) 因此给定 N 和 K，通过式( 18) 和式( 19) 可以很容 易求出相应的低通滤波器． 需要强调说明的是，利用 式( 18) 和式( 19) 构造的分数阶低通滤波器是具有线 性相位特性的( 对称性) ． 2. 2 高通滤波器的构造 构造了具有线性相位特性的低通滤波器之后，将 求解相应的高通滤波器，从而完成图 2 中滤波器组的 构造． 首先，将低通滤波器 H( z) 用多相位方法表示为 H( z) = H0 ( z 2 ) + z － 3H1 ( z 2 ) ，则低通滤波器部分被分解 成两个低通通道 H0 ( z) 和 H1 ( z) ． 为方便使用完全可 重构条件( 2) ～ ( 4) ，再将 H0 ( z) 展开成多相位形式， 即 H0 ( z) = H00 ( z 3 ) + z － 1H01 ( z 3 ) + z － 2H02 ( z 3 ) ，H1 ( z) 和 Gi ( z) ( i = 0，1，2) 按照同样的方式分解，得到完全可重 构条件的等价多相位矩阵表达式: HT ( z － 1 ) GT ( z － 1 [ ) ] H( z) G( z [ ] ) = I3 ． ( 20) 式中 H( z) = H00 ( z) H01 ( z) H02 ( z) H10 ( z) H11 ( z) H12 [ ( z) ]， G( z) = G00 ( z) G01 ( z) G02 ( z) G10 ( z) G11 ( z) G12 ( z) G20 ( z) G21 ( z) G22 ( z        )  ． 通过式( 20) 得到 GT ( z － 1 ) G( z) = I3 － HT ( z － 1 ) H( z) ． ( 21) 下面将通过矩阵谱因式分解的方式求解上式中的 G( z) ． 对于矩阵谱因式分解，文献［13］给出了多种方 法的详细介绍，但这些方法仅限于简单例子的手工计 算． 文献［2］中使用的是对称因式分解方法，但计算结 果还需要最后通过仿酉矩阵调整才能使高通滤波器组 满足近似平移不变特性，因此计算复杂度非常高． 本 文将采用拓普利兹矩阵分解法对式( 21) 进行直接分 解，该方法是一种迭代计算方法，计算方法非常简单． 通过程序计算对比发现，本文采用拓普利兹矩阵 分解法求高通滤波器组具有以下三点优势: ( 1) 对称 因式分解法［2］和使用拓普利兹矩阵分解法迭代 5 000 次计算时间开销虽然接近，但是本文方法在经过 3 000 次拓普利兹矩阵分解迭代运算( 时间为 5 000 次迭代 的 1 /2) 之后，其计算结果误差范围已经在 0. 005 5% ～ 0. 009 7% ． ( 2) 利用拓普利兹矩阵分解法得到的高通 滤波器组已经满足近似平移不变特性，无需再如文献 ［2］一样做额外的计算，因此省去繁杂的调整过程，降 低了计算复杂度． ( 3) 根据不同的数据精度需求和滤 波器的设计复杂程度，采用拓普利兹矩阵分解法也具 有更好的灵活度，同时也更方便计算机实现． 表 1 是文献［2］N = 3 和 K = 1 高通滤波器组． 对 比表 2 中由本文中提出的分解法经过 3 000 次迭代运 算之后的结果，可以看出两种方法求出的高通滤波器 组是近似一样的． 这验证了本文方法的正确性． 表 1 文献［2］中的 N = 3 和 K = 1 的高通滤波器系数 Table 1 High-pass filter coefficients of N = 3 and K = 1 in Ｒef．［2］ g0 ( n) g1 ( n) g2 ( n) 0. 649 177 785 057 41 0 0 － 0. 482 626 543 662 26 0. 637 708 687 474 35 0 － 0. 150 591 301 199 69 － 0. 466 878 032 128 12 0. 645 206 315 833 16 － 0. 014 771 352 175 28 － 0. 148 151 753 049 68 － 0. 490 989 226 274 25 － 0. 001 188 588 020 16 － 0. 022 678 902 296 56 － 0. 131 494 909 897 32 0 0 － 0. 021 526 275 468 89 0 0 － 0. 001 195 904 192 72 · 183 ·