基于加权递推滤波算法的高超声速飞行器鲁棒控制方法

1.本发明属于飞行器制导与控制技术领域，具体涉及一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法。


背景技术：

2.航天器再入过程中经过制导系统优化指令后，会将制导指令传输给控制系统。控制系统在设计过程中，主要涉及到传感器、比较器、执行机构、被控对象、测量反馈等几部分组成。由控制系统将得到的姿态角指令传输到姿态控制器。姿态控制器将收到的控制指令与飞行器实际姿态通过比较器进行对比分析，在控制器计算得到姿态角速率控制指令。姿态角速率控制器再与飞行器状态对比分析，得到控制力矩指令。将控制力矩分配给舵面或者rcs执行机构，来实现对机体姿态的控制。而传感器会将执行之后飞机的状态再次传输给控制器，以此来进行反馈控制。
3.目前大多控制方法的设计只是简单加入一定的干扰，未考虑到大气风场对控制系统的干扰情况。风场变化会对飞行器的执行机构产生一定的干扰，影响飞行器姿态角速率的变化，从而影响姿态角的变化，使其与预期的控制值产生一定的偏差。所以针对风场干扰问题，除了要解决制导律的设计，还需要考虑如何改进控制方法，使其能够更好地适应风场变化，提升整体性能。


技术实现要素：

4.发明目的：为提升控制系统的鲁棒性以及姿态角控制跟踪性能，本发明提供一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法。
5.技术方案：本发明提供一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，包括如下步骤：
6.(1)基于高超声速飞行器获取再入目标任务；
7.(2)通过飞行器自身的制导算法，得到符合目标任务的制导指令；
8.(3)姿态角控制器接收到制导指令之后，通过内部控制器进行计算，结合当前飞行器的姿态角的实时状态量以及干扰观测器和滤波器得到的干扰，对比分析得出满足姿态角控制指令的最新的姿态角速率控制信号；
9.(4)姿态角速率控制器得到由姿态角控制器给出的姿态角速率控制信号后，通过内部控制器进行计算，结合当前飞行器的姿态角速率的实时状态量以及干扰观测器得到的干扰；对比分析得出满足姿态角速率控制指令的最新的控制力矩；
10.(5)在舵面执行器收到控制力矩信号后，根据信号分配至飞行器舵面执行器，飞行器舵面调整控制飞行器姿态；
11.(6)风场干扰后，高超声速飞行器的实时姿态角以及实时姿态角速率会发生变化，反馈系统将受风场扰动的姿态角以及姿态角速率信号传递给干扰观测器，干扰观测器通过将理想控制量与实际状态量对比，得到扰动误差量，并传递给加权递推改进滤波器；
12.(7)基于加权递推改进滤波算法对风紊流产生的控制偏差进行估计，得到精确估计后的姿态角速率干扰以及姿态角干扰，将风场干扰传递给姿态角控制器与姿态角速率控制器，形成闭环。
13.进一步地，步骤(1)所述再入目标任务为飞行器再入初始条件、末端约束、过程约束。
14.进一步地，步骤(3)所述姿态角控制器计算方式为：
[0015][0016]
xc＝[αc,βc,σc]
t
、xs＝[α,β,σ]
t
[0017]
其中，ωc为控制器算法得到的姿态角速率理想控制值；gs为姿态角控制回路系统矩阵；fs为关于姿态角控制回路状态向量的非线性函数；ks＝diag{k
si
}n×n为姿态角回路设计参数，k
si
为正值，为姿态角控制回路常值风场干扰的估计项；xc为控制姿态角，xs为实时姿态角；为控制姿态角的微分项，为控制姿态角估计值，为实时姿态角估计值；αc,βc,σc为迎角、侧滑角以及倾侧角的理想控制值，α,β,σ为反馈系统给出的迎角、侧滑角以及倾侧角的实际值。
[0018]
进一步地，其特征在于，步骤(4)所述姿态角速率控制器计算公式为：
[0019][0020]
x
c(f)
＝[pc,qc,rc]
t
、xf＝[p,q,r]
t
[0021]
其中，mc为控制器算法得到的力矩理想控制值；gf表示姿态角速率控制回路系统矩阵；ff表示关于姿态角速率控制回路状态向量的非线性函数；ωc为姿态角速率的控制值，为姿态角速率的理想控制值的微分项；kf＝diag{k
fi
}n×n为姿态角速率回路设计参数，k
fi
为正值，为姿态角速率控制回路常值风场干扰的估计项；x
c(f)
为控制姿态角速率，xf为实时姿态角速率；为实时姿态角速率的估计值；pc,qc,rc为滚转角速率、俯仰角速率以及偏航角速率的理想控制值，p,q,r为反馈系统给出的滚转角速率、俯仰角速率以及偏航角速率的实际值。
[0022]
进一步地，步骤(6)所述实时姿态角以及实时姿态角速率为：
[0023][0024][0025]
其中，为实时姿态角的微分项，为实时姿态角速率的微分项；fs表示关于姿态角控制回路状态向量的非线性函数，gs表示姿态角控制回路系统矩阵，ff表示关于姿态角速率控制回路状态向量的非线性函数，gf表示姿态角速率控制回路系统矩阵，mc＝[lc,mc,nc]
t
表示控制力矩向量；ds和df分别为姿态角控制回路受到的风场干扰和姿态角速率控制回路受到的风场干扰。
[0026]
进一步地，步骤(6)所述姿态角与姿态角速率干扰观测器的计算处理方法为：
[0027]
[0028]
进一步地，所述步骤(7)实现过程如下：
[0029]
加权递推滤波器的加权递推算法是将递推样本长度为h的样本权重分配按照线性进行分布，斜率为k，区间内第n-i个样本权重a(n-i)的具体表达式为：
[0030][0031]
将该样本值乘以对应权重，得到新的样本值zi表达式为：
[0032]
zi＝si×a(n-i)
[0033]
其中表示姿态角以及姿态角速率干扰观测器得出的风场扰动估计值；
[0034]
滤波算法在容积卡尔曼滤波滤波的基础上进行算法优化：
[0035]
改进前观测量为zi，改进后的观测量为假定置信区间内干扰呈线性等差分布，计算方式为：
[0036][0037]
式中表示自z
k-h
至z
k-1
之间所有观测量求和后的平均值，表达式为：
[0038][0039]
式中
△zk-1|k-2
表示z在k-1与k-2时刻在置信区间内的均值之差；
△zk-1|k-2
表达式为：
[0040][0041]
根据前一时刻z
k-1
的估计值和方差p
k-1|k-1
得到容积点：
[0042][0043]
共有六个容积点，容积点传播公式为：
[0044]
χ
i,k|k-1
＝f(χ
i,k-1|k-1
)
[0045]
依据预测和预测方差p
k|k-1
，计算容积点：
[0046][0047]
根据的计算方程，进行容积点传播：
[0048]zi,k|k-1
＝h(χ
i,k|k-1
),
[0049]
计算量测预测z
k|k-1
、量测预测误差方差sk、状态与量测互协方差ck、增益kk；将测量变量zk按概率分布成与zk/h之和，得到滤波输出以及方差为：
[0050]
[0051]
得到的代替控制器算法中的与
[0052]
有益效果：与现有技术相比，本发明的有益效果：在容积卡尔曼滤波的基础上，进一步引入加权递推算法，调整滤波过程中各个采样点的权重情况，设计非线性干扰观测器与控制器；本发明整合出最佳的权重系数与置信区间和测量误差之间的关系表，为飞行器再入控制系统的误差观测提供更好的权重系数方案，提升控制系统的鲁棒性以及姿态角控制跟踪性能。
附图说明
[0053]
图1是本发明的流程图；
[0054]
图2是最佳权重系数k的分布情况图；
[0055]
图3是采用本发明后控制系统对姿态角速率的跟踪情况效果图；
[0056]
图4是采用本发明后控制系统对姿态角的跟踪情况效果图；
[0057]
图5是采用本发明进行60组仿真后得到的迎角误差均方值；
[0058]
图6是采用本发明进行60组仿真后得到的倾侧角误差均方值。
具体实施方式
[0059]
下面结合附图对本发明做进一步详细说明。
[0060]
本发明提供一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，如图1所示，具体包括以下步骤：
[0061]
步骤1：基于高超声速飞行器获取再入目标任务，包括一般为飞行器再入初始条件、末端约束、过程约束等任务指标。
[0062]
步骤2：通过飞行器自身的制导算法，得到符合目标任务的制导指令αc、βc、σc。
[0063]
步骤3：姿态角控制器接收到制导指令αc、βc、σc之后，通过内部的控制器进行计算，结合当前飞行器的姿态角的实时状态量α、β、σ以及干扰观测器和滤波器得到的干扰对比分析得出满足姿态角控制指令的最新的姿态角速率控制信号pc、qc、rc。
[0064]
姿态角控制器计算方式为：
[0065][0066]
xc＝[αc,βc,σc]
t
、xs＝[α,β,σ]
t
[0067]
其中，ωc为控制器算法得到的姿态角速率理想控制值；gs为姿态角控制回路系统矩阵；fs为关于姿态角控制回路状态向量的非线性函数；ks＝diag{k
si
}n×n为姿态角回路设计参数，k
si
为正值，为姿态角控制回路常值风场干扰的估计项；xc为控制姿态角，xs为实时姿态角；为控制姿态角的微分项，为控制姿态角估计值，为实时姿态角估计值；αc,βc,σc为迎角、侧滑角以及倾侧角的理想控制值，α,β,σ为反馈系统给出的迎角、侧滑角以及倾侧角的实际值。变量上面加“·”为微分项，变量上面加“^”为估计值项。
[0068]
步骤4：姿态角速率控制器得到由姿态角控制器给出的角速率控制指令pc、qc、rc后，通过内部控制器进行计算，结合当前飞行器的姿态角速率的实时状态量p、q、r以及干扰观测器得到的干扰对比分析独处满足姿态角速率控制指令的最新的控制力矩lc、mc、nc。
[0069][0070]
x
c(f)
＝[pc,qc,rc]
t
、xf＝[p,q,r]
t
[0071]
其中，mc为控制器算法得到的力矩理想控制值；gf表示姿态角速率控制回路系统矩阵；ff表示关于姿态角速率控制回路状态向量的非线性函数；ωc为姿态角速率的控制值，为姿态角速率的理想控制值的微分项；kf＝diag{k
fi
}n×n为姿态角速率回路设计参数，k
fi
为正值，为姿态角速率控制回路常值风场干扰的估计项；x
c(f)
为控制姿态角速率，xf为实时姿态角速率；为实时姿态角速率的估计值；pc,qc,rc为滚转角速率、俯仰角速率以及偏航角速率的理想控制值，p,q,r为反馈系统给出的滚转角速率、俯仰角速率以及偏航角速率的实际值。
[0072]
步骤5：在舵面执行器收到控制力矩信号后，根据信号分配至飞行器舵面执行器，飞行器舵面调整控制飞行器姿态。
[0073]
步骤6：风场干扰后，高超声速飞行器的实时姿态角以及实时姿态角速率会发生变化，反馈系统将受风场扰动的姿态角以及姿态角速率信号传递给干扰观测器，干扰观测器通过将理想控制量与实际状态量对比，得到扰动误差量ss与sf，并传递给加权递推改进滤波器。
[0074]
实时姿态角以及实时姿态角速率为：
[0075][0076][0077]
其中，为实时姿态角的微分项，为实时姿态角速率的微分项；fs表示关于姿态角控制回路状态向量的非线性函数，gs表示姿态角控制回路系统矩阵，ff表示关于姿态角速率控制回路状态向量的非线性函数，gf表示姿态角速率控制回路系统矩阵，mc＝[lc,mc,nc]
t
表示控制力矩向量；ds和df分别为姿态角控制回路受到的风场干扰和姿态角速率控制回路受到的风场干扰。
[0078]
姿态角与姿态角速率干扰观测器的计算处理方法为：
[0079]ss
＝x
s-x
c(s) sf＝x
f-x
c(f)
[0080][0081][0082]
步骤7：基于加权递推改进滤波算法对风紊流产生的控制偏差进行估计，得到精确估计后的姿态角速率干扰以及姿态角干扰，将风场干扰传递给姿态角控制器与姿态角速率控制器，形成闭环。
[0083]
加权递推滤波器的加权递推算法是将递推样本长度为h的样本权重分配按照线性进行分布，斜率为k，区间内第n-i个样本权重a(n-i)的具体表达式为：
[0084][0085]
将该样本值乘以对应权重，得到新的样本值zi表达式为：
[0086]
zi＝si×a(n-i)
[0087]
其中表示姿态角以及姿态角速率干扰观测器得出的风场扰动估计值；
[0088]
滤波算法在容积卡尔曼滤波滤波的基础上进行算法优化：
[0089]
改进前观测量为zi，改进后的观测量为假定置信区间内干扰呈线性等差分布，计算方式为：
[0090][0091]
式中表示自z
k-h
至z
k-1
之间所有观测量求和后的平均值，表达式为：
[0092][0093]
式中
△zk-1|k-2
表示z在k-1与k-2时刻在置信区间内的均值之差。
△zk-1|k-2
表达式为：
[0094][0095]
根据前一时刻z
k-1
的估计值和方差p
k-1|k-1
得到容积点：
[0096][0097]
根据步骤(6)中ds和df算法，由于该系统状态为3维状态方程，而容积点个数为维度的2倍，所以共有6个容积点，容积点传播公式为：
[0098]
χ
i,k|k-1
＝f(χ
i,k-1|k-1
),i＝1,2,3，4,5,6
[0099]
计算预测和预测方差：
[0100][0101]
根据此刻的和方差p
k|k-1
计算容积点：
[0102][0103]
根据的计算方程，进行容积点传播：
[0104]zi,k|k-1
＝h(χ
i,k|k-1
),i＝1,2,3，4,5,6
[0105]
进行量测预测z
k|k-1
、量测预测误差方差(新息方差)sk、状态与量测互协方差ck、增益kk的表达式为：
[0106][0107]
将测量变量zk按概率分布成与zk/h之和。得到滤波输出以及方差为：
[0108][0109]
以此滤波方法得到的代替控制器算法中的与
[0110]
本发明主要目的是研究风干扰下高超声速再入飞行器的姿态控制系统的设计问题，因此观测跟踪风场干扰带来的姿态角误差是实现优化的前提条件因。本发明采用状态量与控制量跟踪误差的均方值来评价跟踪性能的好坏。考虑到风场具备随机性、不同飞行器的采样周期不一、置信区间不一，进行重复实验取最小误差均方值对应的权重系数k作为最佳值。不同置信区间以及采样误差情况下的仿真得到的最佳权重系数k值如表1：
[0111]
表1最佳权重分布表
[0112]
[0113][0114]
如图2所示，为蒙特卡洛仿真试验后，与表1相对应的更直观展现最佳系数k设置图，能够看出在不同置信区间以及采样误差下，取最小跟踪误差时对应权重系数k的值。
[0115]
如图3、图4所示，为采用本发明的单次实验得到的再入期间姿态角速率、姿态角的
控制值与常规控制方法受风场干扰的实际值之间对比效果，能够看出改进后的跟踪误差相较于改进之前，误差幅值更小。姿态角速率、姿态角跟踪性能更好。
[0116]
如图5、图6所示，为经过60组蒙特卡洛打靶仿真后，本方法改进后的实验得到的再入期间姿态角速率、姿态角的跟踪误差均方值。能够看出改进后的姿态角速率跟踪误差均方值稳定在0.05以内，改进后的姿态角跟踪误差均方值稳定在0.3以内。
[0117]
此外，应当理解，虽然本说明书按照实施例加以描述，但并非每个实施例仅包含一个独立的技术方案，说明书的这种叙述方式仅仅是为清楚起见，本领域技术人员应当将说明书作为一个整体，各实施例中的技术方案也可以经适当组合，形成本领域技术人员可以理解的其他实施方式。

技术特征：
1.一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，包括以下步骤：(1)基于高超声速飞行器获取再入目标任务；(2)通过飞行器自身的制导算法，得到符合目标任务的制导指令；(3)姿态角控制器接收到制导指令之后，通过内部控制器进行计算，结合当前飞行器的姿态角的实时状态量以及干扰观测器和滤波器得到的干扰，对比分析得出满足姿态角控制指令的最新的姿态角速率控制信号；(4)姿态角速率控制器得到由姿态角控制器给出的姿态角速率控制信号后，通过内部控制器进行计算，结合当前飞行器的姿态角速率的实时状态量以及干扰观测器得到的干扰；对比分析得出满足姿态角速率控制指令的最新的控制力矩；(5)在舵面执行器收到控制力矩信号后，根据信号分配至飞行器舵面执行器，飞行器舵面调整控制飞行器姿态；(6)风场干扰后，高超声速飞行器的实时姿态角以及实时姿态角速率会发生变化，反馈系统将受风场扰动的姿态角以及姿态角速率信号传递给干扰观测器，干扰观测器通过将理想控制量与实际状态量对比，得到扰动误差量，并传递给加权递推改进滤波器；(7)基于加权递推改进滤波算法对风紊流产生的控制偏差进行估计，得到精确估计后的姿态角速率干扰以及姿态角干扰，将风场干扰传递给姿态角控制器与姿态角速率控制器，形成闭环。2.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，步骤(1)所述再入目标任务为飞行器再入初始条件、末端约束、过程约束。3.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，步骤(3)所述姿态角控制器计算方式为：x
c
＝[α
c
,β
c
,σ
c
]
t
、x
s
＝[α,β,σ]
t
其中，ω
c
为控制器算法得到的姿态角速率理想控制值；g
s
为姿态角控制回路系统矩阵；f
s
为关于姿态角控制回路状态向量的非线性函数；k
s
＝diag{k
si
}
n
×
n
为姿态角回路设计参数，k
si
为正值，为姿态角控制回路常值风场干扰的估计项；x
c
为控制姿态角，x
s
为实时姿态角；为控制姿态角的微分项，为控制姿态角估计值，为实时姿态角估计值；α
c
,β
c
,σ
c
为迎角、侧滑角以及倾侧角的理想控制值，α,β,σ为反馈系统给出的迎角、侧滑角以及倾侧角的实际值。4.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，其特征在于，步骤(4)所述姿态角速率控制器计算公式为：x
c(f)
＝[p
c
,q
c
,r
c
]
t
、x
f
＝[p,q,r]
t
其中，m
c
为控制器算法得到的力矩理想控制值；g
f
表示姿态角速率控制回路系统矩阵；f
f
表示关于姿态角速率控制回路状态向量的非线性函数；ω
c
为姿态角速率的控制值，为姿态角速率的理想控制值的微分项；k
f
＝diag{k
fi
}
n
×
n
为姿态角速率回路设计参数，k
fi
为正
值，为姿态角速率控制回路常值风场干扰的估计项；x
c(f)
为控制姿态角速率，x
f
为实时姿态角速率；为实时姿态角速率的估计值；p
c
,q
c
,r
c
为滚转角速率、俯仰角速率以及偏航角速率的理想控制值，p,q,r为反馈系统给出的滚转角速率、俯仰角速率以及偏航角速率的实际值。5.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，步骤(6)所述实时姿态角以及实时姿态角速率为：特征在于，步骤(6)所述实时姿态角以及实时姿态角速率为：其中，为实时姿态角的微分项，为实时姿态角速率的微分项；f
s
表示关于姿态角控制回路状态向量的非线性函数，g
s
表示姿态角控制回路系统矩阵，f
f
表示关于姿态角速率控制回路状态向量的非线性函数，g
f
表示姿态角速率控制回路系统矩阵，m
c
＝[l
c
,m
c
,n
c
]
t
表示控制力矩向量；d
s
和d
f
分别为姿态角控制回路受到的风场干扰和姿态角速率控制回路受到的风场干扰。6.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，步骤(6)所述姿态角与姿态角速率干扰观测器的计算处理方法为：s
s
＝x
s-x
c(s)
s
f
＝x
f-x
c(f)c(f)
7.根据权利要求1所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，所述步骤(7)实现过程如下：加权递推滤波器的加权递推算法是将递推样本长度为h的样本权重分配按照线性进行分布，斜率为k，区间内第n-i个样本权重a(n-i)的具体表达式为：将该样本值乘以对应权重，得到新的样本值z
i
表达式为：z
i
＝s
i
×
a
(n-i)
其中表示姿态角以及姿态角速率干扰观测器得出的风场扰动估计值；滤波算法在容积卡尔曼滤波滤波的基础上进行算法优化：改进前观测量为z
i
，改进后的观测量为假定置信区间内干扰呈线性等差分布，计算方式为：式中表示自z
k-h
至z
k-1
之间所有观测量求和后的平均值，表达式为：
式中
△
z
k-1|k-2
表示z在k-1与k-2时刻在置信区间内的均值之差；
△
z
k-1|k-2
表达式为：根据前一时刻z
k-1
的估计值和方差p
k-1|k-1
得到容积点：容积点传播公式为：χ
i,k|k-1
＝f(χ
i,k-1|k-1
)依据预测和预测方差p
k|k-1
，计算容积点：根据的计算方程，进行容积点传播：z
i,k|k-1
＝h(χ
i,k|k-1
),计算量测预测z
k|k-1
、量测预测误差方差s
k
、状态与量测互协方差c
k
、增益k
k
；将测量变量z
k
按概率分布成与z
k
/h之和，得到滤波输出以及方差为：得到的代替控制器算法中的与8.根据权利要求7所述的基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，其特征在于，所述容积点为六个。

技术总结
本发明公开了一种基于加权递推滤波算法的高超声速飞行器鲁棒控制方法，为了精确估算风场对再入飞行器的影响并进行校正，将环境干扰分为稳定干扰与非稳定干扰；在改进容积卡尔曼滤波的基础上，引入了加权递推算法，优化非线性干扰观测器与控制器；通过加权递推方法，调整滤波过程中各个采样点的权重；根据仿真结果统计整合出最佳的权重系数与置信区间和测量误差之间的关系表，为飞行器再入控制系统的误差观测提供更好的权重系数方案。本发明提升了滤波实时性，通过非线性干扰观测器与鲁棒控制器实现降低跟踪误差的均方值，提高了飞行器再入指令跟踪能力。再入指令跟踪能力。再入指令跟踪能力。


技术研发人员：唐明明 都延丽 刘燕斌 魏昀鹏 王文凯
受保护的技术使用者：南京航空航天大学
技术研发日：2023.04.06
技术公布日：2023/7/7