技术领域
[0001] 本
发明属于控制工程与科学的技术领域，具体涉及一种自适应最优相角陷波滤波器设计方法。
背景技术
[0002] 在实际所用到的伺服机构中，往往广泛存在着一些集中在窄频带的干扰
信号。国内外很多专家学者也对此提出了很多的滤波或者控制方法，比如在低频段(100–600Hz)，常常采用峰值滤波器来有效地抑制低频窄带干扰，然而，峰值滤波器却不能抑制在开环交叉
频率周围的中频干扰，因为它的固有
相位会损耗相位裕度，在扰动频率附近会使得系统回路灵敏度增益失真。而在中频段(1.6kHz左右)，
现有技术常常采用一个相位超前峰值滤波器来抑制系统中出现的中频段窄带干扰，这种滤波器为了保证相位裕度，增加了一个额外的π/2的相位超前的改进微分器将灵敏度曲线平滑化。对于
硬盘伺服系统中的高频(4-10kHz)窄带干扰抑制，多采用
锁相伺服
控制器来抑制由悬架振动引起的偏差干扰。这些滤波器或控制器设计虽然有很好的抑制窄带干扰的能
力，但是它们都只能在有限的
频率范围有效地抑制窄带干扰。自适应滤波器的出现大大促进了未知频率干扰领域的发展，但是一般的自适应滤波器方法中存在的主要问题是：
算法复杂，滤波参数维数过高，自适应时间太长，系统实时性大大降低。
发明内容
[0003] 本发明针对多个频率变化的未知单频
干扰信号和伺服系统本身存在的不确定性等随机干扰在内的多源干扰信号对系统模型的影响，为了实现高
精度控制，提供一种自适应最优相角陷波滤波器设计方法。
[0004] 本发明提供的一种自适应最优相角陷波滤波器设计方法，应用系统回路模
块包括：对象模型模块，鲁棒控制器模块和自适应最优相角陷波滤波器模块。鲁棒控制器模块和对象模型模块
串联在系统回路的输入端与输出端之间，自适应最优相角陷波滤波器模块并联在系统回路的输入端与鲁棒控制器模块之间，同时串联
开关K。
[0005] 信号传递过程：一方面自适应最优相角陷波滤波器模块的
输入信号x1(n)与误差信号e(n)一起进入自适应最优相角陷波滤波器模块中进行自适应滤波，输出自适应最优相角陷波滤波器模块的
输出信号y1(n)，另一方面x1(n)信号与y1(n)信号进行
叠加得到鲁棒控制器模块的输入信号x2(n)，x2(n)信号经过鲁棒控制器模块和对象模型模块完成控制过程，得到对象模型模块的输出信号y2(n)，y2(n)信号与单频干扰信号d(n)叠加，得到误差信号e(n)，误差信号e(n)一方面与x1(n)信号进入自适应最优相角陷波滤波器模块中，另一方面叠加宽频随机干扰信号，被输入信号r(n)减去，进行
负反馈，得到自适应最优相角陷波滤波器模块的输入信号x1(n)，构成整体的
闭环系统。
[0006] 本发明提供的一种自适应最优相角陷波滤波器设计方法，具体实施步骤包括：
[0007] 步骤1、通过频率响应实验或物理公式推导建立对象模型模块，用传递函数表示为P(s)，离散化后表示为P(z)。
[0008] 步骤2、鲁棒控制器模块的传递函数表达式C(s)，离散化后表示为C(z)的设计方法如下：
[0009] ①设原系统支路灵敏度函数
[0010] 补灵敏度函数T0(s)＝1-S0(s)
[0011] ②性能加权函数W1(s)，表示成如下形式：
[0012]
[0013] K1为高增益，-ε为二重极点，取值范围为0<ε<1；
[0014] ③鲁棒加权函数W2(s)
[0015]
[0016] 其中Δ代表对象模型幅值相对误差，j为虚数单位，ωi'为频率，i=1，…，m，m表示实验重复次数，n表示每次实验的
采样点数，k表示第k次实验，(Mik,φik)为频率ωi下系统响应的幅值和相角，(Mi,φk)为该频率点对应的拟合后的P(s)中的幅值和相角。
[0017] W2(s)的传递函数表达式设计为二阶模型
[0018]
[0019] 其中K2为鲁棒权值函数增益，ζ1,ζ2分别为W2(s)零点多项式和极点多项式的阻尼比，ζ1,ζ2∈(0,1)，ω1,ω2分别为W2(s)零点多项式和极点多项式的固有频率；
[0020] ④得到W1(s)、W2(s)、S0(s)和T0(s)之后，调用矩阵实验室
软件Matrix Laboratory中鲁棒控制器工具箱中的混合灵敏度控制器函数，指令为mixsyn，寻找正则的鲁棒控制器，使得W1(s)、W2(s)、S0(s)和T0(s)满足无穷范数形式的权值条件：
[0021]
[0022] 其中γ是H∞的性能指标，取值范围是γ＞0，正则的含义是鲁棒控制器传递函数分母中s的次数要高于分子中s的次数；
[0023] 得到鲁棒控制器连入闭环系统回路中，既能保证系统的鲁棒
稳定性，又能达到抑制宽频域的随机干扰信号的效果。
[0024] 步骤3、自适应最优相角陷波滤波器模块的传递函数表达式F(s)，离散化后表示为F(z)的设计方法如下：
[0025]
[0026] 其中ω0为需要抑制干扰的频率，s＝jω0时相角 ζ是阻尼比ζ∈(0,1)，Kf是最优相角滤波器正增益。
[0027] 步骤4、得到z变换之后对象模型模块的传递函数表达式P(z)，鲁棒控制器模块的传递函数表达式C(z)和自适应最优相角陷波滤波器模块的传递函数表达式F(z)以后，系统回路的信号过程如下：
[0028] （1）自适应最优相角陷波滤波器模块的输入信号x1(n)一方面与误差信号e(n)一起进入自适应最优相角陷波滤波器模块中，进行自适应滤波，直到输出自适应最优相角陷波滤波器模块的输出信号y1(n)；具体步骤如下：
[0029] ①对于未知频率的单频干扰信号，需要采用根据对象模型改进的最小均方值（LMS）算法，进行部分滤波器参数的自适应过程：
[0030]
[0031] c0,c1,......,cM为C(z)·P(z)的分子系数，构成定常矩阵C的元素；M为自然数，d1,......,dN为C(z)·P(z)的分母系数，构成定常矩阵D的元素；N为正整数。
[0032] ②F(s)离散化之后，自适应最优相角陷波滤波器离散化传递函数F(z)在n时刻的表达式为：
[0033]
[0034] 其中b0，b1，b2为F(z)的分子系数，构成定常矩阵B的元素；a1(n)，a2为F(z)的分母系数，构成矩阵A(n)的元素；其中A(n)与n时刻有关，a1(n)随n时刻变化而变化。
[0035] ③进行自适应滤波后，此时从z0到z-M算起，共M+1个系数分别构成了C(z)·P(z)和F(z)的分子和分母系数矩阵，在n时刻表示为如下形式：
[0036]
[0037] 自适应最优相角陷波滤波器模块的输入信号x1(n)和鲁棒控制器模块的输入信号x2(n)从n时刻往前计算到n-M时刻的值分别为状态矩阵X1(n)、X2(n)，自适应最优相角陷波滤波器模块的输出信号y1(n)和对象模型模块的输出信号y2(n)从n-1时刻往前计算到n-M-1时刻的值分别为状态矩阵Y1(n)、Y2(n)，有如下关系式：
[0038] Xi(n)＝[xi(n)…xi(n-M)]1×(M+1),i＝1,2
[0039] Yi(n)＝[yi(n-1)…yi(n-M-1)]1×(M+1),i＝1,2
[0040] ④由F(z)的表达式得到n时刻的离散状态形式：
[0041]
[0042] 进而得到：
[0043] (1+a1(n)z-1+a2z-2)y1(n)＝(b0+b1z-1+b2z-2)x1(n)
[0044] 又由z-kx(n)＝x(n-k)，z-ky(n)＝y(n-k)得
[0045] y1(n)+a1(n)y1(n-1)+a2y1(n-2)＝b0x1(n)+b1x1(n-1)+b2x1(n-2)
[0046] 矩阵B与矩阵X1(n)的转置相结合后，得到结果为
[0047] BX1T(n)=b0x1(n)+b1x1(n-1)+b2x1(n-2)
[0048] 矩阵A(n)与矩阵Y1(n)的转置相结合后，得到结果为
[0049] A(n)Y1T(n)=a1(n)y1(n-1)+a2y1(n-2)
[0050] BX1T(n)的值与A(n)Y1T(n)的值进行相减操作得到y1(n)
[0051] y1(n)＝BX1T(n)-A(n)Y1T(n)
[0052] 其中T代表矩阵的转置。
[0053] 即
[0054] y1(n)＝b0x1(n)+b1x1(n-1)+b2x1(n-2)-a1(n)y1(n-1)-a2y1(n-2)
[0055] （2）x1(n)信号另一方面与y1(n)信号进行叠加得到鲁棒控制器模块的输入信号x2(n)，n时刻x1(n),x2(n),y1(n)数值关系
[0056] x2(n)＝x1(n)+y1(n)
[0057] 从而状态矩阵X1(n),X2(n),Y1(n)有如下关系：
[0058] X2(n)＝X1(n)+Y1(n+1)
[0059] （3）x2(n)信号经过鲁棒控制器模块和对象模型模块完成控制过程得到对象模型模块的输出信号y2(n)，具体操作步骤同得到自适应最优相角陷波滤波器模块的输出信号y1(n)的步骤相同。
[0060] 矩阵C与矩阵X2(n)的转置相结合后，得到结果为
[0061] CX2T(n)=CX1T(n)+CY1T(n+1)
[0062] 选取ck为矩阵C中按z降幂排序中第一个不为零的系数，0≤k＜N
[0063] 得到：
[0064] CX1T(n)=ckx1(n-k)+…+cMx1(n-M)
[0065] CY1T(n+1)＝cky1(n-k)+…+cMy1(n-M)
[0066] 矩阵D与矩阵Y2(n)的转置相结合后，得到结果为：
[0067] DY2T(n)＝d1y2(n-1)-…-dNy2(n-N-1)
[0068] 由C(z)·P(z)表达式得到离散状态形式
[0069] y2(n)＝CX2T(n)-DY2T(n)
[0070] y2(n)=ckx1(n-k)+…+cMx1(n-M)+cky1(n-k)+…+cMy1(n-M)-d1y2(n-1)-…-dNy2(n-N-1)
[0071] （4）对象模型模块的输出信号y2(n)与单频干扰信号d(n)叠加，得到误差信号e(n)为：
[0072] e(n)＝d(n)+y2(n)
[0073] （5）误差信号e(n)一方面与自适应最优相角陷波滤波器模块的输入信号x1(n)进入自适应最优相角陷波滤波器模块中，另一方面叠加宽频随机干扰信号，被输入信号r(n)减去，进行负反馈，构成整体的闭环系统。具体步骤如下：
[0074] ①取误差信号e(n)的均方值为代价函数J(n)：
[0075] J(n)＝E[e2(n)]
[0076] ②代价函数J(n)对于自适应参数a1(n-k)的梯度为：
[0077]
[0078] 表示梯度，a1(n-k)为自适应参数；N为D矩阵中非零元素的个数。
[0079] ③按照最小均方值算法，采用瞬时值估计梯度矢量：
[0080]
[0081] ^表示瞬时值；
[0082] ④采用改进型LMS
迭代算法，求出 信号对自适应参数 的梯度：
[0083]
[0084] ⑤由偏导数性质，先求出 的偏导数，得到运算结果：
[0085]
[0086] ⑥由偏导数性质，再求出 的偏导数，得到运算结果：
[0087]
[0088] ⑦最小均方值算法的瞬时值：
[0089]
[0090] ⑧则自适应参数迭代公式为：
[0091]
[0092] μ为迭代收敛步长，满足 λmax是自适应最优相角陷波滤波器模块的输出信号y1(n)构成的
自相关矩阵R的最大特征值。
[0093] ⑨将 的表达式代入上式得到自适应参数迭代公式的最终表达形式：
[0094]
[0095] 将n+1时刻自适应最优相角陷波滤波器系数a1(n+1)取 的值继续运算，到某一个时间ΔT之后，ΔT的大小与收敛步长μ有关，μ越小ΔT越大，达到抑制未知单频干扰信号的效果，自适应最优相角陷波滤波器设计完成。
[0096] 本发明的优点和积极效果在于：
[0097] 1)本发明提出了一种自适应最优相角陷波滤波器的设计方法，其滤波器的
零点能根据补灵敏度函数在干扰频率处的相角变化而改变，能够在干扰频率处形成陷波，消除已知频率和未知频率的单频干扰信号，设计方法简单有效。
[0098] 2)本发明研究了一种控制器和滤波器的并联实现结构，进而引入H∞混合灵敏度控制设计鲁棒控制器，从而实现了对宽频随机干扰信号的最优抑制；
[0099] 3)本发明基于特殊的最优相角陷波滤波器的模型及其与系统的并联方式，创造性地将其与最小均方值自适应滤波算法结合起来，使陷波滤波器的陷波中心频率能
跟踪未知单频干扰信号频率，达到对未知单频干扰信号的抑制作用；
[0100] 4)本发明在设计最优相角自适应陷波滤波器的过程中，巧妙发现模型中的特点，简化自适应参数为一个，大大缩短自适应时间，增强了实时性，并在一定精度范围内满足抗干扰条件。
[0101] 5）本发明提出的一种自适应最优相角陷波滤波器设计方法，它能最低限度地影响闭环系统的稳定性，而且在无限频率范围内通过分配零阶滤波器来有效地抑制未知单频干扰信号，并配合设计H∞混合灵敏度鲁棒控制器，以提高系统的鲁棒性及随机干扰抑制性能，使系统兼具随机干扰抑制能力与未知频率单频干扰信号抑制能力。
附图说明
[0102] 图1为本发明自适应最优相角陷波滤波器设计方法的
流程图。
[0103] 图2为本发明的自适应最优相角陷波滤波器的数字实现流程图。
具体实施方式
[0104] 下面将结合附图和具体
实施例对本发明作进一步地详细说明。
[0105] 本发明提出的一种自适应最优相角陷波滤波器设计方法，可以应用于具有宽频随机干扰信号和多个未知频率单频干扰信号等多干扰源信号同时存在的情况下，达到最优抑制。
[0106] 所述一种自适应最优相角陷波滤波器设计方法所应用的系统回路，如图1所示，所应用的系统回路模块包括：对象模型模块，鲁棒控制器模块以及自适应最优相角陷波滤波器模块。鲁棒控制器模块和对象模型模块串联在系统回路的输入端与输出端之间。自适应最优相角陷波滤波器模块并联在输入端与鲁棒控制器模块之间，同时串联开关K。
[0107] 鲁棒控制器模块的功能是控制闭环系统鲁棒稳定并且抑制伺服系统本身存在的不确定性相关的随机干扰信号，对象模型模块的功能是用传递函数的方法表述被控对象，自适应最优相角陷波滤波器模块的功能是对未知单频干扰信号的抑制。
[0108] 首先，在对象数学模型传递函数P(s)
基础之上，通过调整性能加权函数W1(s)和鲁棒加权函数W2(s)，得出合适的鲁棒控制器传递函数C(s)；其次，设计自适应最优相角陷波滤波器模块的传递函数F(s)，自适应最优相角陷波滤波器模块并联接入部分系统回路中；最后，根据P(s)·C(s)的传递函数，针对滤波器参数a1采用改进型最小均方值算法，使陷波中心频率跟踪未知单频干扰信号频率。
[0109] 信号传递方面如图1所示；一方面自适应最优相角陷波滤波器模块的输入信号x1(n)与误差信号e(n)一起进入自适应最优相角陷波滤波器模块中进行自适应滤波，输出自适应最优相角陷波滤波器模块的输出信号y1(n)，另一方面x1(n)信号与y1(n)信号进行叠加得到鲁棒控制器模块的输入信号x2(n)，x2(n)信号经过鲁棒控制器模块和对象模型模块完成控制过程，得到对象模型模块的输出信号y2(n)，y2(n)信号与单频干扰信号d(n)叠加，得到误差信号e(n)，误差信号e(n)一方面与x1(n)信号进入自适应最优相角陷波滤波器模块中，另一方面叠加宽频随机干扰信号，被输入信号r(n)减去，进行负反馈，构成整体的闭环系统。
[0110] 本发明提供的一种自适应最优相角陷波滤波器设计方法，具体实施步骤包括：
[0111] 步骤1、通过频率响应实验或物理公式推导建立对象的数学模型，也就是对象模型模块，拉普拉斯变换后对象模型模块的传递函数表达式为P(s)，传递函数P(s)z变换之后得到离散传递函数表达式为P(z)；其中s为拉普拉斯变换算子。
[0112] 鲁棒控制器模块的传递函数表达式为C(s)，自适应最优相角陷波滤波器模块的传递函数表达式为F(s)，z变换之后鲁棒控制器模块的传递函数表达式为C(z)，自适应最优相角陷波滤波器模块的传递函数表达式为F(z)，z代表z变换的算子。
[0113] 步骤2、鲁棒控制器模块拉普拉斯变换后的传递函数表达式C(s)，z变换离散化后表示为C(z)的设计方法如下：
[0114] ①在得到对象模型模块之后，运用相关现代控制理论如H2和H∞最优控制技术，使控制对象鲁棒稳定，以H∞混合灵敏度鲁棒控制器设计为例，对于对象和控制器回路，设原系统支路灵敏度函数为：
[0115]
[0116] 补灵敏度函数为：
[0117]
[0118] ②W1(s)为性能加权函数，反映了原系统灵敏度函数在低频段的性能特征及形状要求，从而在低频段获得高增益，增强抗干扰的抑制能力，W1(s)选择如下形式：
[0119]
[0120] K1为高增益，-ε代表二重极点，取0<ε<1；
[0121] ③W2(s)为鲁棒加权函数，因为鲁棒加权函数是不确定函数，由高频段动态和建模参数不确定性所决定，反映了平台系统本身的性能。
[0122] 如果对象模型模块是稳定的，并且它的传递函数是通过频率响应的实验得到的，测定在一系列的频率下的对象响应的幅值和相角，这一系列频率为ωi'，i=1，…，m，m表示实验重复次数，另外，实验重复n次，得到的对象模型传递函数P(s)是扫频之后图像拟合而成，扫频的每一个频率点对应的幅值和相角设为(Mik,φik)，同样输入情况下该频率点对应的拟合后的P(s)图像中的幅值和相角设为(Mi,φk)，W2(s)设计的条件满足：
[0123]
[0124] j代表虚数单位，Δ代表对象模型幅值相对误差，n表示每次实验的采样点数，k表示第k次实验。
[0125] 通过这个方法设计出来的W2(s)能让|W2(jωi')|在所需频率段内将Δ表达式的图线完全包络，其中W2(s)的传递函数表达式设计为二阶模型：
[0126]
[0127] K2为鲁棒权值函数增益，ζ1,ζ2分别为W2(s)零点多项式和极点多项式的阻尼比，ζ1,ζ2∈(0,1)，ω1,ω2分别为W2(s)零点多项式和极点多项式的固有频率，ζ1ω1s,ζ2ω2s分别代表三个系数乘积；
[0128] ④得到W1(s)、W2(s)、S0(s)和T0(s)之后，调用Matrix Laboratory(矩阵实验室软件)中鲁棒控制器工具箱中的混合灵敏度控制器函数，指令为mixsyn，混合灵敏度问题是寻找正则的鲁棒控制器，使得W1(s)、W2(s)、S0(s)和T0(s)足无穷范数形式的权值条件：
[0129]
[0130] 其中γ是H∞的性能指标，取值范围是γ＞0，正则的含义是鲁棒控制器传递函数分母中s的次数要高于分子中s的次数；
[0131] ⑤鲁棒控制器模块进行拉普拉斯变换后的传递函数表达式C(s)，z变换离散化之后为C(z)，这样在对数幅频特性曲线中，设计出的原系统支路灵敏度函数S0(s)在低于交叉频率的频段内幅值小于0dB（分贝），就能达到抑制宽频域的随机干扰信号的效果；
[0132] 步骤3、自适应最优相角陷波滤波器模块拉普拉斯变换后的传递函数表达式F(s)，z变换离散化后表示为F(z)的设计方法如下：
[0133] （1）自适应最优相角陷波滤波器模块的传递函数表达式F(s)，离散化后表示为F(z)的表达式：
[0134]
[0135] 其中ω0为需要抑制干扰的频率，ζ是阻尼比ζ∈(0,1)，Kf是最优相角滤波器正增益，当s＝jω0时相角
[0136] （2）自适应最优相角陷波滤波器模块抑制已知单频干扰信号的过程如下：
[0137] ①自适应最优相角陷波滤波器模块通过并联结构连入系统中，得出系统总体灵敏度函数S(s)表达式：
[0138] 分子分母同时乘上（1+P(s)C(s)）可得：
[0139]
[0140] 打开开关K，从输入到输出的传递函数表示为原系统支路灵敏度函数S0(s)[0141]
[0142] 合上开关K，滤波器支路从输入到输出的传递函数表示为滤波器支路灵敏度函数SF(s)，具体地，
[0143]
[0144] 最后系统总体灵敏度函数S(s)串行解耦成原系统支路灵敏度函数S0(s)与滤波器支路灵敏度函数SF(s)的乘积：
[0145] S(s)＝S0(s)·SF(s) （11）
[0146] 将干扰频率s＝jω0代入到滤波器支路灵敏度函数SF(s)中得到相关增益比较关系：
[0147]
[0148] |SF(jω0)|的值由|T0(jω0)|和|F(jω0)|决定，由于T0(jω0)的值与鲁棒控制器模块C(s)和对象模型模块P(s)有关，已确定，所以只需考虑F(jω0)的值。
[0149] ②自适应最优相角陷波滤波器模块的传递函数表达式F(s)，离散化后表示为F(z)的表达式：
[0150]
[0151] 其中Kf是最优相角陷波滤波器的正增益，ω0为需要抑制干扰的频率，ζ是阻尼比ζ∈(0,1)， 与补灵敏度函数T0(s)有关，s＝jω0时相角 两个实零点分别为z1＝0，
[0152] 当相角 满足如下关系时，滤波器支路灵敏度函数的幅值|SF(jω0)|达到最小值：
[0153]
[0154] 由
[0155] Re(F(jω0))表示对函数F(jω0)取复数的
实部，Im(F(jω0))表示对函数F(jω0)取复数的
虚部，
[0156] 可知相角：
[0157]
[0158] 所以
[0159]
[0160]
[0161] Kf与ζ的取值在s＝jω0的相角满足arg[T0(jω0)]+arg[F(jω0)]＝0时，|SF(jω0)|达到最小值，|SF(jω0)|min在对数幅频特性图像上远远低于0dB，进而滤波器支路灵敏度函数SF(s)与原系统支路灵敏度函数S0(s)相乘之后，就能在单频干扰信号频率处对原系统支路灵敏度函数S0(s)产生一个比较大的衰减，并且在其他频率处几乎不改变原系统支路灵敏度函数S0(s)的形状，而得到最后的系统总体灵敏度函数S(s)，达到抑制已知单频干扰信号的效果。
[0162] 对于已知频率的单频干扰信号，设计相应的最优相角陷波滤波器来实现抑制作用，而对于未知频率的单频干扰信号，则需要采用根据对象模型模块改进的最小均方值（LMS）算法，对最优相角滤波器进行部分滤波器参数的自适应滤波过程，使得陷波中心频率能跟踪未知单频干扰信号频率，达到对未知频率单频干扰信号的抑制作用。
[0163] 步骤4、得到z变换之后对象模型模块的传递函数表达式P(z)，鲁棒控制器模块的传递函数表达式C(z)和自适应最优相角陷波滤波器模块的传递函数表达式F(z)以后，系统回路的信号过程如下：
[0164] （1）如图1所示：自适应最优相角陷波滤波器模块的输入信号x1(n)一方面与误差信号e(n)一起进入自适应最优相角陷波滤波器模块中，进行自适应滤波，直到输出自适应最优相角陷波滤波器模块的输出信号y1(n)，自适应最优相角陷波滤波器的数字实现流程步骤如下：
[0165] ①改进型最小均方值算法在鲁棒控制器模块及对象模型模块中的实现。为简化算法，考虑鲁棒控制器模块离散化传递函数和对象模型模块离散化传递函数之积：
[0166]
[0167] 其中c0,c1,......,cM为鲁棒控制器离散化传递函数和对象模型离散化传递函数之积C(z)·P(z)的分子系数，该分子系数构成定常矩阵C的元素；M为自然数，d1,......,dN为C(z)·P(z)的分母系数，该分母系数构成定常矩阵D的元素；N为正整数。z0,......z-M共M+1个系数为矩阵C，D的列。
[0168] ②确定自适应参数个数。将自适应最优相角陷波滤波器传递函数离散化之后得到：
[0169]
[0170] 其中，b0，b1，b2代表离散化传递函数F(z)的分子系数，a1，a2代表离散化传递函数F(z)的分母系数。
[0171] 通过实验方法可以验证自适应最优相角陷波滤波器的一个特点：该滤波器零点的变化虽然对滤波器灵敏度函数的形状有一定的影响，但都是在陷波中心频率两边的一些频率段有相应灵敏度增益的增加或者衰减，对于陷波中心频率改变不大，对陷波中心频率处的灵敏度增益也几乎没有改变，即对最优相角陷波滤波器抗干扰性能特征影响不大。所以在一定精度允许范围内，自适应参数可以不考虑滤波器零点参数的相关变化，也就是滤波器分子中的系数b0，b1，b2可以在某一频率认定为不变的，而对于滤波器分母系数a1，a2来说，对于只考虑干扰频率ω0改变的情况下，a2不变化，也就是说，对于自适应最优相角滤波器算法来说，只需要考虑a1一个参数的变化。
[0172] 自适应最优相角陷波滤波器传递函数离散化之后，如图2所示，自适应最优相角陷波滤波器离散化传递函数F(z)在n时刻的表达式为：
[0173]
[0174] b0，b1，b2构成定常矩阵B的元素；a1(n)，a2为F(z)的分母系数，构成矩阵A(n)的元素；其中A(n)与n时刻有关，a1(n)随n时刻变化而变化。a1(0)，a1(1)，......,a1(k)设为初始值；当n-k时刻得到的自适应参数设为a1(n-k)；z0，……，z-2为矩阵A(n)，B的列。
[0175] ③进行自适应滤波后，从z0到z-M算起，共M+1个系数分别构成了C(z)·P(z)和F(z)的分子和分母系数矩阵，其中A(n)与n时刻有关，B、C、D均为定常矩阵，在n时刻表示为如下形式：
[0176]
[0177] 自适应最优相角陷波滤波器模块的输入信号x1(n)和鲁棒控制器模块的输入信号x2(n)从n时刻往前计算到n-M时刻的值分别为状态矩阵分别为X1(n)和X2(n)，自适应最优相角陷波滤波器模块的输出信号y1(n)和对象模型模块的输出信号y2(n)从n-1时刻往前计算到n-M-1时刻的值分别为状态矩阵Y1(n)和Y2(n)，有如下关系式：
[0178] Xi(n)＝[xi(n)…xi(n-M)]1×(M+1)(i＝1,2) (21)
[0179] Yi(n)＝[yi(n-1)…yi(n-M-1)]1×(M+1)(i＝1,2) (22)
[0180] ④由F(z)的表达式可得n时刻的离散状态形式：
[0181]
[0182] 即：
[0183] (1+a1(n)z-1+a2z-2)y1(n)＝(b0+b1z-1+b2z-2)x1(n)
[0184] 又由z-kx(n)＝x(n-k)，z-ky(n)＝y(n-k)
[0185] 可得：
[0186] y1(n)+a1(n)y1(n-1)+a2y1(n-2)＝b0x1(n)+b1x1(n-1)+b2x1(n-2) （24）[0187] 如图2所示，矩阵B与矩阵X1(n)的转置相结合后，得到结果为
[0188] BX1T(n)=b0x1(n)+b1x1(n-1)+b2x1(n-2) （25）
[0189] 矩阵A(n)与矩阵Y1(n)的转置相结合后，得到结果为
[0190] A(n)Y1T(n)=a1(n)y1(n-1)+a2y1(n-2) （26）
[0191] ⑤BX1T(n)的值与A(n)Y1T(n)的值进行相减操作得到y1(n)
[0192] y1(n)＝BX1T(n)-A(n)Y1T(n) (27)
[0193] 其中T代表矩阵的转置。
[0194] 即
[0195] y1(n)＝b0x1(n)+b1x1(n-1)+b2x1(n-2)-a1(n)y1(n-1)-a2y1(n-2) （28）[0196] （2）自适应最优相角陷波滤波器模块的输入信号x1(n)另一方面与自适应最优相角陷波滤波器模块的输出信号y1(n)进行叠加得到鲁棒控制器模块的输入信号x2(n)，[0197] x1(n),x2(n)和y1(n)数值关系，如图1所示：
[0198] x2(n)＝x1(n)+y1(n) （29）
[0199] 从而状态矩阵X1(n),X2(n),Y1(n)有如下关系：
[0200] X2(n)＝X1(n)+Y1(n+1)
[0201] 进行矩阵转置后得到：
[0202] X2T(n)=X1T(n)+Y1T(n+1) （30）
[0203] （3）x2(n)信号经过鲁棒控制器模块和对象模型模块完成控制过程得到对象模型模块的输出信号y2(n)，得到对象模型模块的输出信号y2(n)的具体操作步骤，同得到自适应最优相角陷波滤波器模块的输出信号y1(n)的步骤相同。
[0204] 矩阵C与矩阵X2(n)的转置相结合后，得到结果为：
[0205] CX2T(n)=CX1T(n)+CY1T(n+1) （31）
[0206] 选取ck为矩阵C中按z降幂排序中第一个不为零的系数，0≤k＜N
[0207] 得到：
[0208] CX1T(n)=ckx1(n-k)+…+cMx1(n-M) （32）
[0209] CY1T(n+1)＝cky1(n-k)+…+cMy1(n-M) （33）
[0210] 所以：
[0211] CX2T(n)=ckx1(n-k)+…+cMx1(n-M)+cky1(n-k)+…+cMy1(n-M) （34）[0212] 矩阵D与矩阵Y2(n)的转置相结合后，得到结果为：
[0213] DY2T(n)＝d1y2(n-1)-…-dNy2(n-N-1) （35）
[0214] 由C(z)·P(z)表达式可得到离散状态形式：
[0215] y2(n)=CX2T(n)-DY2T(n) （36）
[0216] y2(n)＝ckx1(n-k)+…+cMx1(n-M)+cky1(n-k)+…+cMy1(n-M)-d1y2(n-1)-…-dNy2(n-N-1（) 37）
[0217] （4）对象模型模块的输出信号y2(n)与单频干扰信号d(n)叠加，得到误差信号e(n)，如图1所示：n时刻的误差信号e(n)可表示为：
[0218]
[0219] （5）误差信号e(n)一方面与自适应最优相角陷波滤波器模块的输入信号x1(n)进入自适应最优相角陷波滤波器模块中，进行自适应滤波过程。另一方面叠加宽频随机干扰，被输入信号r(n)减去，进行负反馈，构成整体的闭环系统。自适应滤波过程步骤如下：
[0220] ①设误差信号e(n)的均方值代价函数为J(n)：
[0221] J(n)＝E[e2(n)] （39）
[0222] ②代价函数J(n)对于自适应参数a1(n-k)的梯度为：
[0223]
[0224] 该改进算法的目的即满足代价函数J(n)趋向最小值，采用改进型LMS迭代算法，在常规的梯度迭代算法中一般要考虑延时最少相关的误差信号与迭代参数，即n时刻的误差信号对n-q时刻的迭代参数求偏导，用 表示梯度，其中q越小收敛效果越好，设a1(n-k)为自适
[0225] 应参数；N为D矩阵中非零元素的个数。
[0226] ③按照最小均方值算法，采用瞬时值估计梯度矢量的方法：
[0227]
[0228] ^表示瞬时值；
[0229] ④ 信号对自适应参数 的梯度计算过程如下：
[0230] 对误差信号e(n)，CY1T(n+1)与自适应参数a1(n-k)相关，又由于CX1T(n)与a1(n-k)无关，DY2T(n)与a1(n-k-1)有关，即e(n)求偏导时只用考虑CY1T(n+1)对a1(n-k)求偏导，由于CY1T(n+1)中与a1(n-k)有关的项只有cky1(n-k)，即：cky1(n-k)对y1(n-k)的梯度计算过程：
[0231] 由n时刻的e(n)的展开式，与a1(n-k)有关的项只有y1(n-k)，故由二阶偏微分性质可得：
[0232]
[0233] ⑤由偏导数性质，先求出 对 的偏导数，得到运算结果
[0234]
[0235]
[0236] ⑥自适应最优相角陷波滤波器模块的输出信号y1(n)，经过k+1次延时环节之后，得到y1(n-k-1)，
[0237] 如附图2所示，结合公式z-ka1(n)＝a1(n-k)可得：
[0238] z-(k+1)y1(n)＝y1(n-k-1)
[0239] 由偏导数性质，再求出 对 的偏导数，得到运算结果：
[0240]
[0241] ⑦最小均方值算法的瞬时值：
[0242]
[0243] ⑧如附图2所示，得到y1(n)后，设定前k+1个时刻的值a1(0)，a1(1)，......,a1(k)为初始值，从第k+1时刻开始，每隔k+1进行一次迭代计算，当第n-k时刻的自适应参数a1(n-k)进行迭代计算后，对应的自适应参数值为a1(n+1)；
[0244] 自适应参数a1(n-k)与k+1倍的单次参数修正值相加，得到n+1时刻的自适应参数值a1(n+1)，其中单次参数修正值为迭代收敛步长μ乘以二分之一倍的负的代价函数的梯度，则自适应参数迭代公式为：
[0245]
[0246] 其中，μ为迭代收敛步长，满足 λmax是输出信号y1(n)构成的自相关矩阵R的最大特征值。
[0247] ⑨将 的表达式代入上式可得自适应参数迭代公式的最终表达形式：
[0248]
[0249] 将n+1时刻自适应最优相角陷波滤波器系数a1(n+1)取 的值继续运算，到某一个时间ΔT之后，ΔT的大小与收敛步长μ有关，μ越小，ΔT越大，自适应最优相角陷波滤波器的陷波中心频率会跟踪上单频干扰信号频率，自适应过程完成，从而达到抑制单频干扰信号的效果，自适应最优相角陷波滤波器设计完成。
[0250] 该设计方法重点研究了一种特殊的最优相角陷波滤波器的结构和自适应算法，并以此为基础研究了一种控制器和滤波器的并联实现结构，进而引入H∞混合灵敏度控制设计鲁棒控制器，从而实现了对随机干扰信号和未知单频干扰信号信号的最优抑制。这一算法可以应用于具有宽频随机干扰，多个未知频率的单频干扰信号等多源干扰同时存在情况的最优抑制。
[0251] 本发明
说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。
