CN102735259B - 一种基于多层状态估计器的卫星控制系统故障诊断方法 - Google Patents

Abstract

一种基于多层状态估计器的卫星控制系统故障诊断方法，能够在硬件冗余条件不满足的情况下，实现光学敏感器和惯性敏感器的故障隔离。该方法基本流程为：首先，建立用于描述光学敏感器和惯性敏感器之间输出关系的解析模型；然后，基于模型分别设计两层状态估计器，第1层状态估计器产生的残差对光学敏感器故障和惯性敏感器故障都敏感，而第2层状态估计器的残差不受惯性敏感器故障的影响。这样，通过检测第1层状态估计器的残差并与事先设定的阈值相比较，可以判断是否发生故障；如果发生故障，可以根据第2层状态估计器的残差区分光学敏感器故障和惯性敏感器故障。

Description

一种基于多层状态估计器的卫星控制系统故障诊断方法

技术领域

本发明属于卫星姿态控制领域，涉及一种基于多层状态估计器的卫星控制系统故障诊断方法。

背景技术

现代航天器大多结构复杂，工作环境恶劣，运行时间长，在轨工作过程中极易发生故障。一旦发生故障，若不采取积极有效的措施进行应对，很可能导致航天任务失败，从而带来巨大的资源浪费和不利的社会影响。

作为构成卫星的各个分系统中最关键也是最复杂的一个分系统，卫星姿态控制系统的可靠性是星上多个分系统正常运行的基本保障。如遥测遥控分系统要实现与地面站的通信，要求天线指向正确的位置；能源分系统要保障正常的电力供应，要求太阳帆板跟踪太阳；对于对地观测卫星而言，有效载荷的正常工作要求卫星进行高精度高稳定度姿态控制。因此，卫星姿态控制系统能否正常工作往往关系到整个航天任务的成败。

卫星姿态控制系统由姿态敏感器、控制器、执行机构和卫星本体一起构成闭环控制回路。基于姿态敏感器测量和确定卫星相对于空间某些已知基准目标的方位；控制器对测得的信息进一步处理后确定卫星姿态，并根据所确定的姿态按满足设计要求的控制率发出指令，控制执行机构按控制指令产生所需的控制力矩，实现卫星姿态控制。卫星姿态控制系统的故障率较高且危害较大。

为了提高卫星姿态控制系统的可靠性，一方面，应提高其组成部件如各敏感器、执行机构和控制器等自身的可靠性，降低发生故障的可能性；另一方面，应对关键部件采用冗余设计，并通过故障诊断技术消除或减弱故障造成的不利影响。

卫星姿态敏感器包括光学敏感器和惯性敏感器，其中，光学敏感器有太阳敏感器、地球敏感器等，惯性敏感器以陀螺为代表。一般采用一致性检验的方法进行卫星姿态敏感器故障诊断，如通过比较两个光学敏感器输出的方式进行故障检测，但当硬件冗余条件不满足时(如无故障光学敏感器个数少于2个)，不能采用这种方法进行故障诊断。为了解决硬件冗余条件不满足情况下卫星姿态敏感器的故障诊断问题，本发明提出基于多层状态估计器的光学敏感器和惯性敏感器故障诊断方法。

发明内容

本发明的技术解决问题是：克服现有技术的不足，提供了一种基于多层状态估计器的卫星姿态敏感器故障诊断方法。能够在硬件冗余条件不满足的情况下，实现光学敏感器和惯性敏感器的故障隔离。

本发明的技术解决方案是：

一种基于多层状态估计器的卫星控制系统故障诊断方法，步骤如下：

(1)建立用于描述星载的光学敏感器和惯性敏感器之间输出关系的故障诊断模型，即卫星俯仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型，

卫星俯仰通道故障诊断模型形式如下所示：

x_p(k+1)＝A_px_p(k)+B_pu_p(k)

y_p(k)＝C_px_p(k)

其中，k表示离散的时间，x_p(k)表示第k步的俯仰通道状态量，A_p表示俯仰通道状态矩阵，B_p表示俯仰通道输入矩阵，u_p(k)表示俯仰通道输入量，y_p(k)表示俯仰通道观测量，C_p表示俯仰通道观测矩阵；

卫星滚动-偏航通道故障诊断模型形式如下所示：

x_r(k+1)＝A_rx_r(k)+B_ru_r(k)

y_r(k)＝C_rx_r(k)

其中，x_r(k)表示第k步的滚动-偏航通道状态量，A_r表示滚动-偏航通道状态矩阵，B_r表示滚动-偏航通道输入矩阵，u_r(k)表示滚动-偏航通道输入量，y_r(k)表示滚动-偏航通道观测量，C_r示滚动-偏航通道观测矩阵；

(2)根据步骤(1)中的俯卫星仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型设计基于俯仰通道故障诊断模型的第1层状态估计器和基于滚动-偏航通道故障诊断模型的第1层状态估计器；

基于俯仰通道故障诊断模型的第1层状态估计器形式为：

${\hat{x}}_p(k+1|k)=A_p{\hat{x}}_p\left(k\right|k)+B_p{u}_p\left(k\right)$

${\hat{x}}_p(k+1|k+1)=$

${\hat{x}}_p(k+1|k)+K_p(k+1)[y_p(k+1)-C_p{\hat{x}}_p(k+1|k)]$

其中，和分别表示俯仰通道状态量的预测值和估计值，K_p(k+1)表示俯仰通道的滤波增益阵；

根据俯仰通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_p\left(k\right)\right|\left|\right)$

其中，N表示故障检测窗口的长度，r_p(k)表示俯仰通道第1层状态估计器的测量残差，其形式如下所示：

$r_p\left(k\right)=y_p\left(k\right)-C_p{\hat{x}}_p\left(k\right|k-1)$

基于滚动-偏航通道故障诊断模型的第1层状态估计器为：

${\hat{x}}_r(k+1|k)=A_r{\hat{x}}_r\left(k\right|k)+B_r{u}_r\left(k\right)$

${\hat{x}}_r(k+1|k+1)=$

${\hat{x}}_r(k+1|k)+K_r(k+1)[y_r(k+1)-C_r{\hat{x}}_r(k+1|k)]$

其中，和分别表示滚动-偏航通道状态量的预测值和估计值，K_r(k+1)表示滚动-偏航通道的滤波增益阵；

根据滚动-偏航通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_r\left(k\right)\right|\left|\right)$

其中，r_p(k)表示滚动-偏航通道第1层状态估计器的测量残差，其形式如下所示：

$r_r\left(k\right)=y_r\left(k\right)-C_r{\hat{x}}_r\left(k\right|k-1)$

(3)如果LD_p(k)≥T_D或LD_r(k)≥T_D，则卫星姿态控制系统中的光学敏感器或惯性敏感器发生了故障，进入步骤(4)；反之则光学敏感器和惯性敏感器均未发生故障；其中，T_D为预设的故障检测阈值；

(4)设计第2层状态估计器；

基于俯仰通道故障诊断模型的第2层状态估计器形式为：

z_p(k+1)＝F_pz_p(k)+G_py_p(k)+T_pB_pu_p(k)

e_p(k)＝K_pz_p(K)+P_py_p(k)

其中，z_p(k)表示俯仰通道第2层状态估计器的状态量，e_p(k)表示俯仰通道第2层状态估计器的输出，F_p、G_p、T_p、K_p和P_p表示俯仰通道第2层状态估计器设计参数矩阵，其中，

F_p＝ρ_pI

I表示单位阵，ρ_p为可调参数且|ρ_p|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_p{A}_p-F_p{T}_p=G_p{C}_p\\ K_p{T}_p+P_p{C}_p=0\\ K_p{T}_p{B}_p=0\end{array}\right.$

根据俯仰通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_k\left(k\right)$

基于滚动-偏航通道故障诊断模型设计的第2层状态估计器形式为：

z_r(k+1)＝F_rz_r(k)+G_ry_r(k)+T_rB_ru_r(k)

e_r(k)＝K_rz_r(k)+P_ry_r(k)

其中，z_r(k)表示滚动-偏航通道第2层状态估计器的状态量，e_r(k)表示滚动-偏航通道第2层状态估计器的输出，F_r、G_r、T_r、K_r和P_r表示滚动-偏航通道第2层状态估计器设计参数矩阵，其中，

F_r＝ρ_rI

ρ_r为可调参数且|ρ_r|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_r{A}_r-F_r{T}_r=G_r{C}_r\\ K_r{T}_r+P_r{C}_r=0\\ K_r{T}_r{B}_r=0\end{array}\right.$

根据滚动-偏航通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_r\left(k\right)$

(5)如果LI_p(k)≥T_I或LI_r(k)≥T_I，则光学敏感器发生故障，反之，则惯性敏感器发生故障，其中，T_I表示预设的故障隔离阈值；

所述俯仰通道状态量x_p(k)的形式为：

$x_p\left(k\right)=\left[\begin{array}{c}\mathrm{\θ}\left(k\right)\\ d_p\left(k\right)\\ b_p\left(k\right)\\ b_{\mathrm{\θ}}\left(k\right)\end{array}\right]$

其中，θ(k)表示卫星俯仰角，d_p(k)表示俯仰轴陀螺指数漂移误差，b_p(k)表示俯仰轴陀螺常值漂移误差，b_θ(k)表示地球敏感器俯仰角观测量常值偏差；

所述俯仰通道状态矩阵A_p的形式为：

$A_p=\left[\begin{array}{cccc}1& -1& -1& 0\\ 0& 1-{\mathrm{\α}}_p& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]$

其中，α_p表示俯仰轴陀螺指数偏移误差系数。

所述俯仰通道输入矩阵B_p的形式为：

$B_p=\left[\begin{array}{c}1\\ 0\\ 0\\ 0\end{array}\right].$

所述俯仰通道输入量u_p(k)的形式为：

u_p(k)＝ω_o(k)+g_p(k)

其中，ω_o(k)表示卫星轨道角速率，g_p(k)表示俯仰轴陀螺测量输出；

所述俯仰通道观测量形式为：

$y_p\left(k\right)=\left[\begin{array}{c}{\mathrm{\θ}}_h\left(k\right)\\ {\mathrm{\θ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，θ_h(k)和θ_sh(k)分别表示俯仰轴地球敏感器和太阳敏感器测量输出；

所述俯仰通道观测矩阵形式为： $C_p=\left[\begin{array}{cccc}1& 0& 0& 1\\ 1& 0& 0& 0\end{array}\right].$

所述滚动-偏航通道状态量的形式为：

$x_r\left(k\right)=\left[\begin{array}{c}\mathrm{\φ}\left(k\right)\\ \mathrm{\ψ}\left(k\right)\\ d_r\left(k\right)\\ d_y\left(k\right)\\ b_r\left(k\right)\\ b_y\left(k\right)\\ b_{\mathrm{\φ}}\left(k\right)\end{array}\right]$

其中，φ(k)和ψ(k)分别表示卫星滚动角和偏航角，d_r(k)和d_y(k)分别表示滚动轴和偏航轴陀螺指数漂移误差，b_r(k)和b_y(k)分别表示滚动轴和偏航轴陀螺常值漂移误差，b_φ(k)表示地球敏感器滚动角观测量常值偏差；

所述滚动-偏航通道状态矩阵形式为：

$A_r=\left[\begin{array}{ccccccc}1& {\mathrm{\ω}}_0& -1& 0& -1& 0& 0\\ -{\mathrm{\ω}}_0& 1& 0& -1& 0& -1& 0\\ 0& 0& 1-{\mathrm{\α}}_r& 0& 0& 0& 0\\ 0& 0& 0& 1-{\mathrm{\α}}_y& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 1\end{array}\right]$

其中，α_r和α_y分别表示滚动轴和偏航轴陀螺指数偏移误差系数；

所述滚动-偏航通道输入矩阵形式为：

$B_r=\left[\begin{array}{cc}1& 0\\ 0& 1\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\end{array}\right]$

所述滚动-偏航通道输入量形式为：

$u_r\left(k\right)=\left[\begin{array}{c}{g}_r\left(k\right)\\ g_y\left(k\right)\end{array}\right]$

其中，g_r(k)和g_y(k)分别表示滚动轴和偏航轴陀螺测量输出；

所述滚动-偏航通道观测量形式为：

$y_r\left(k\right)=\left[\begin{array}{c}{\mathrm{\φ}}_h\left(k\right)\\ {\mathrm{\φ}}_{\mathrm{sh}}\left(k\right)\\ {\mathrm{\ψ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，φ_h(k)和φ_sh(k)分别表示滚动轴地球敏感器和太阳敏感器测量输出，ψ_sh(k)表示太阳敏感器滚动轴测量输出；

所述滚动-偏航通道观测矩阵形式为：

$C_r=\left[\begin{array}{ccccccc}1& 0& 0& 0& 0& 0& 1\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\end{array}\right]$

本发明与现有技术相比的有益效果是：

国内卫星多采用基于硬件冗余的卫星姿态敏感器故障诊断方式，但当硬件冗余条件不满足时(如光学敏感器个数少于2个，陀螺个数少于4个)，不能采用这种方法进行故障诊断。针对这一问题，本发明提出一种基于多层状态估计器的卫星姿态敏感器故障诊断方法，该方法能够在硬件冗余条件不满足的情况下，实现光学敏感器和惯性敏感器的故障检测和故障隔离。

附图说明

图1为本发明流程图；

图2a为无故障情况下俯仰滤波器残差曲线；

图2b为无故障情况下俯仰隔离观测器残差曲线；

图2c为无故障情况下滚偏滤波器残差曲线；

图2d为无故障情况下滚偏隔离观测器残差曲线；

图3a为光学敏感器故障情况下俯仰滤波器残差曲线；

图3b为光学敏感器故障情况下俯仰隔离观测器残差曲线

图3c为光学敏感器故障情况下滚偏滤波器残差曲线；

图3d为光学敏感器故障情况下滚偏隔离观测器残差曲线；

具体实施方式

下面结合附图对本发明的具体实施方式进行进一步的详细描述。

如图1所示，本发明提出一种基于多层状态估计器的卫星姿态敏感器故障诊断方法，步骤如下：

(1)建立用于描述光学敏感器和惯性敏感器之间输出关系的故障诊断模型，即卫星俯仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型。假设卫星处于三轴稳定的工作模式下，此时卫星的姿态角和姿态角速度都是小量，卫星的俯仰姿态和滚动-偏航姿态可以解耦。其中，卫星俯仰通道故障诊断模型形式如下所示：

x_p(k+1)＝A_px_p(k)+B_pu_p(k)

y_p(k)＝C_px_p(k)

其中，k表示离散的时间，x_p(k)表示第k步的俯仰通道状态量，A_p表示俯仰通道状态矩阵，B_p表示俯仰通道输入矩阵，u_p(k)表示俯仰通道输入量，y_p(k)表示俯仰通道观测量，C_p表示俯仰通道观测矩阵。俯仰通道状态量的形式为

$x_p\left(k\right)=\left[\begin{array}{c}\mathrm{\θ}\left(k\right)\\ d_p\left(k\right)\\ b_p\left(k\right)\\ b_{\mathrm{\θ}}\left(k\right)\end{array}\right]$

其中，θ(k)表示卫星俯仰角，d_p(k)表示俯仰轴陀螺指数漂移误差，b_p(k)表示俯仰轴陀螺常值漂移误差，b_θ(k)表示地球敏感器俯仰角观测量常值偏差。状态矩阵形式为

$A_p=\left[\begin{array}{cccc}1& -1& -1& 0\\ 0& 1-{\mathrm{\α}}_p& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]$

其中，α_p表示俯仰轴陀螺指数偏移误差系数。输入矩阵形式为

$B_p=\left[\begin{array}{c}1\\ 0\\ 0\\ 0\end{array}\right]$

输入量形式为

u_p(k)＝ω_o(k)+g_p(k)

其中，ω_o(k)表示卫星轨道角速率，g_p(k)表示俯仰轴陀螺测量输出。观测量形式为

$y_p\left(k\right)=\left[\begin{array}{c}{\mathrm{\θ}}_h\left(k\right)\\ {\mathrm{\θ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，θ_h(k)和θ_sh(k)分别表示俯仰轴地球敏感器和太阳敏感器测量输出。观测矩阵形式为

$C_p=\left[\begin{array}{cccc}1& 0& 0& 1\\ 1& 0& 0& 0\end{array}\right]$

滚动-偏航通道故障诊断模型形式如下所示：

x_r(k+1)＝A_rx_r(k)+B_ru_r(k)

y_r(k)＝C_rx_r(k)

其中，x_r(k)表示第k步的滚动-偏航通道状态量，A_r表示滚动-偏航通道状态矩阵，B_r表示滚动-偏航通道输入矩阵，u_r(k)表示滚动-偏航通道输入量，y_r(k)表示滚动-偏航通道观测量，C_r表示滚动-偏航通道观测矩阵。状态量的形式为

$x_r\left(k\right)=\left[\begin{array}{c}\mathrm{\φ}\left(k\right)\\ \mathrm{\ψ}\left(k\right)\\ d_r\left(k\right)\\ d_y\left(k\right)\\ b_r\left(k\right)\\ b_y\left(k\right)\\ b_{\mathrm{\φ}}\left(k\right)\end{array}\right]$

其中，φ(k)和ψ(k)分别表示卫星滚动角和偏航角，d_r(k)和d_r(k)分别表示滚动轴和偏航轴陀螺指数漂移误差，b_r(k)和b_y(k)分别表示滚动轴和偏航轴陀螺常值漂移误差，b_φ(k)表示地球敏感器滚动角观测量常值偏差。状态矩阵形式为

$A_r=\left[\begin{array}{ccccccc}1& {\mathrm{\ω}}_0& -1& 0& -1& 0& 0\\ -{\mathrm{\ω}}_0& 1& 0& -1& 0& -1& 0\\ 0& 0& 1-{\mathrm{\α}}_r& 0& 0& 0& 0\\ 0& 0& 0& 1-{\mathrm{\α}}_y& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 1\end{array}\right]$

其中，α_r和α_y分别表示滚动轴和偏航轴陀螺指数偏移误差系数。输入矩阵形式为

$B_r=\left[\begin{array}{cc}1& 0\\ 0& 1\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\end{array}\right]$

输入量形式为

$u_r\left(k\right)=\left[\begin{array}{c}{g}_r\left(k\right)\\ g_y\left(k\right)\end{array}\right]$

其中，g_r(k)和g_y(k)分别表示滚动轴和偏航轴陀螺测量输出。观测量形式为

$y_r\left(k\right)=\left[\begin{array}{c}{\mathrm{\φ}}_h\left(k\right)\\ {\mathrm{\φ}}_{\mathrm{sh}}\left(k\right)\\ {\mathrm{\ψ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，φ_h(k)和φ_sh(k)分别表示滚动轴地球敏感器和太阳敏感器测量输出，ψ_sh(k)表示太阳敏感器滚动轴测量输出。观测矩阵形式为

$C_r=\left[\begin{array}{ccccccc}1& 0& 0& 0& 0& 0& 1\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\end{array}\right]$

上述故障诊断模型可用于两层状态估计器的设计。

(2)根据步骤(1)中的俯卫星仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型设计基于俯仰通道故障诊断模型的第1层状态估计器和基于滚动-偏航通道故障诊断模型的第1层状态估计器；

基于俯仰通道故障诊断模型设计的第1层状态估计器形式如下所示：

${\hat{x}}_p(k+1|k)=A_p{\hat{x}}_p\left(k\right|k)+B_p{u}_p\left(k\right)$

${\hat{x}}_p(k+1|k+1)=$

${\hat{x}}_p(k+1|k)+K_p(k+1)[y_p(k+1)-C_p{\hat{x}}_p(k+1|k)]$

其中，和分别表示俯仰通道状态量的预测值和估计值，可以通过递推解算获得，K_p(k+1)表示俯仰通道的滤波增益阵。

根据俯仰通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_p\left(k\right)\right|\left|\right)$

其中，N表示故障检测窗口的长度，r_p(k)表示俯仰通道第1层状态估计器的测量残差，其形式如下所示：

$r_p\left(k\right)=y_p\left(k\right)-C_p{\hat{x}}_p\left(k\right|k-1)$

基于滚动-偏航通道故障诊断模型设计的第1层状态估计器形式如下所示：

${\hat{x}}_r(k+1|k)=A_r{\hat{x}}_r\left(k\right|k)+B_r{u}_r\left(k\right)$

${\hat{x}}_r(k+1|k+1)=$

${\hat{x}}_r(k+1|k)+K_r(k+1)[y_r(k+1)-C_r{\hat{x}}_r(k+1|k)]$

其中，和分别表示滚动-偏航通道状态量的预测值和估计值，可以通过递推解算获得，K_r(k+1)表示滚动-偏航通道的滤波增益阵。

根据滚动-偏航通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_r\left(k\right)\right|\left|\right)$

其中，r_p(k)表示滚动-偏航通道第1层状态估计器的测量残差，其形式如下所示：

$r_r\left(k\right)=y_r\left(k\right)-C_r{\hat{x}}_r\left(k\right|k-1)$

(3)监测故障检测残差，当故障检测残差超过某个事先给定的故障检测阈值T_D时，则判断卫星姿态控制系统中的光学敏感器或惯性敏感器发生了故障，之后进入步骤(4)；反之，则判断光学敏感器和惯性敏感器均未发生故障。即：如果LD_p(k)≥T_D或LD_r(k)≥T_D，则卫星姿态控制系统中的光学敏感器或惯性敏感器发生了故障，进入步骤(4)；反之则光学敏感器和惯性敏感器均未发生故障；其中，T_D为预设的故障检测阈值。

(4)为了区分光学敏感器故障和惯性敏感器故障，设计第2层状态估计器。基于俯仰通道故障诊断模型设计的第2层状态估计器形式如下所示：

z_p(k+1)＝F_pz_p(k)+G_py_p(k)+T_pB_pu_p(k)

e_p(k)＝K_pz_p(k)+P_py_p(k)

其中，z_p(k)表示俯仰通道第2层状态估计器的状态量，e_p(k)表示俯仰通道第2层状态估计器的输出，F_p、G_p、T_p、K_p和P_p表示俯仰通道第2层状态估计器设计参数矩阵，其中，

F_p＝ρ_pI

I表示适当维数的单位阵，ρ_p为可调参数，要求|ρ_p|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_p{A}_p-F_p{T}_p=G_p{C}_p\\ K_p{T}_p+P_p{C}_p=0\\ K_p{T}_p{B}_p=0\end{array}\right.$

根据俯仰通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_k\left(k\right)$

基于滚动-偏航通道故障诊断模型设计的第2层状态估计器形式如下所示：

z_r(k+1)＝F_rz_r(k)+G_ry_r(k)+T_rB_ru_r(k)

e_r(k)＝K_rz_r(k)+P_ry_r(k)

其中，z_r(k)表示滚动-偏航通道第2层状态估计器的状态量，e_r(k)表示滚动-偏航通道第2层状态估计器的输出，F_r、G_r、T_r、K_r和P_r表示滚动-偏航通道第2层状态估计器设计参数矩阵，其中，

F_r＝ρ_rI

ρ_r为可调参数，要求|ρ_r|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_r{A}_r-F_r{T}_r=G_r{C}_r\\ K_r{T}_r+P_r{C}_r=0\\ K_r{T}_r{B}_r=0\end{array}\right.$

根据滚动-偏航通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_r\left(k\right)$

(5)在检测出光学敏感器或惯性敏感器发生故障的情况下，监测故障隔离残差，如果故障隔离残差超出事先给定的故障隔离阈值T_I，则判断光学敏感器发生故障，反之，则判断惯性敏感器发生故障，即：如果LI_p(k)≥T_I或LI_r(k)≥T_I，则光学敏感器发生故障，反之，则惯性敏感器发生故障，其中，T_I表示预设的故障隔离阈值。T_D的取值范围为0.1～0.15，T_I的取值范围为0.05～0.07。

下面，以某低轨卫星姿态控制系统为诊断对象，通过仿真实例验证本发明所述方法的有效性。卫星姿态敏感器包括陀螺、地球敏感器和太阳敏感器。故障诊断模型参数如下。

俯仰通道：

$A_p=\left[\begin{array}{cccc}1& -1& -1& 0\\ 0& 0.999& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]$

$B_p=\left[\begin{array}{c}1\\ 0\\ 0\\ 0\end{array}\right]$

$C_p=\left[\begin{array}{cccc}1& 0& 0& 1\\ 1& 0& 0& 0\end{array}\right]$

滚动-偏航通道：

$A_r=\left[\begin{array}{ccccccc}1& 0.001& -1& 0& -1& 0& 0\\ -0.001& 1& 0& -1& 0& -1& 0\\ 0& 0& 0.999& 0& 0& 0& 0\\ 0& 0& 0& 0.999& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 1\end{array}\right]$

$B_r=\left[\begin{array}{cc}1& 0\\ 0& 1\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\end{array}\right]$

$C_r=\left[\begin{array}{ccccccc}1& 0& 0& 0& 0& 0& 1\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\end{array}\right]$

为了实现惯性敏感器和光学敏感器的故障检测和故障隔离，设计两层状态估计器。第1层状态估计器为卡尔曼滤波器，卡尔曼滤波增益阵K_p(k+1)和K_r(k+1)的计算方法可参见中国宇航出版社1998年出版的由屠善澄主编的《卫星姿态动力学与控制》一书。卡尔曼滤波器的递推计算过程可参考西北工业大学出版社1998出版的由秦永元、张洪钺、汪叔华编写的《卡尔曼滤波与组合导航原理》一书。第2层状态估计器为观测器，选择观测器设计参数矩阵F_p＝0.618I₂，F_r＝ρ_rI₂

其中，I₂为二阶单位阵。进而，采用本发明所述方法，通过求解矩阵方程可得其它设计参数矩阵。

俯仰通道：

$G_p=\left[\begin{array}{cc}-1& -1\\ 0& 1\end{array}\right]$

$T_p=\left[\begin{array}{cccc}-5.235& -13.734& -13.705& -2.617\\ 2.617& 6.867& 6.852& 0\end{array}\right]$

K_p＝[1  2]

P_p＝[2.618  -2.618]

滚动-偏航通道：

$G_r=\left[\begin{array}{ccc}-1& 0& -1\\ -1& -1& -2\end{array}\right]$

$T_r=-\left[\begin{array}{ccccccc}2.624& 2.610& 6.885& 6.848& 6.871& 6.834& 2.617\\ 5.249& 5.221& 13.771& 13.697& 13.743& 13.668& 2.617\end{array}\right]$

K_r＝[2-1]

P_r＝[2.618  -2.618  0]

选择故障检测窗口的长度N＝10，阈值T_D＝0.1，T_I＝0.05，对以下2种情况进行仿真研究，分别是：

(a)无故障情况；

(b)光学敏感器故障，故障表现为地球敏感器常值偏差由0.9°增大为1.3°。

分别对以上2种情况进行仿真，两层状态估计器产生的残差曲线如图2(图2a，图2b，图2c，图2d)和图3(图3a，图3b，图3c，图3d)所示。图中，上半部分(图2a，图2c，图3a，图3c)为根据第1层状态估计器的状态估计结果产生的故障检测残差，下半部分(图2b，图2d，图3b，图3d)为根据第2层状态估计器的状态估计结果产生的故障隔离残差。

从如图2(图2a，图2b，图2c，图2d)和图3(图3a，图3b，图3c，图3d)不难看出，无故障时两层状态估计器的残差LD_p(k)、LD_r(k)、LI_p(k)和LI_r(k)均保持在0值附近。而第1层状态估计器的残差LD_p(k)和LD_r(k)超过阈值，表明光学敏感器或惯性敏感器发生了故障；第2层状态估计器的残差LI_p(k)和LI_r(k)超过阈值，表明光学敏感器发生故障。

仿真研究表明，所提方法可以实现对光学敏感器和惯性敏感器故障的检测和隔离。相对于现有的基于硬件冗余的敏感器故障诊断方式，该方法对硬件冗余要求低，易于星上实现。本项专利的主要技术内容可用于各类卫星姿态控制系统故障诊断，具有广阔的应用前景。

本发明说明书中未作详细描述的内容属于本领域专业技术人员的公知技术。

Claims (1)

1.一种基于多层状态估计器的卫星控制系统故障诊断方法，其特征在于步骤如下：

(1)建立用于描述星载的光学敏感器和惯性敏感器之间输出关系的故障诊断模型，即卫星俯仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型，

卫星俯仰通道故障诊断模型形式如下所示：

x_p(k+1)＝A_px_p(k)+B_pu_p(k)

y_p(k)＝C_px_p(k)

其中，k表示离散的时间，x_p(k)表示第k步的俯仰通道状态量，A_p表示俯仰通道状态矩阵，B_p表示俯仰通道输入矩阵，u_p(k)表示俯仰通道输入量，y_p(k)表示俯仰通道观测量，C_p表示俯仰通道观测矩阵；

卫星滚动-偏航通道故障诊断模型形式如下所示：

x_r(k+1)＝A_rx_r(k)+B_ru_r(k)

y_r(k)＝C_rx_r(k)

其中，x_r(k)表示第k步的滚动-偏航通道状态量，A_r表示滚动-偏航通道状态矩阵，B_r表示滚动-偏航通道输入矩阵，u_r(k)表示滚动-偏航通道输入量，y_r(k)表示滚动-偏航通道观测量，C_r表示滚动-偏航通道观测矩阵；

(2)根据步骤(1)中的卫星俯仰通道故障诊断模型和卫星滚动-偏航通道故障诊断模型设计基于俯仰通道故障诊断模型的第1层状态估计器和基于滚动-偏航通道故障诊断模型的第1层状态估计器；

基于俯仰通道故障诊断模型的第1层状态估计器形式为：

${\hat{x}}_p(k+1|k)=A_p{\hat{x}}_p\left(k\right|k)+B_p{u}_p\left(k\right)$

$\begin{array}{c}{\hat{x}}_p(k+1|k+1)=\\ {\hat{x}}_p(k+1|k)+K_p(k+1)[y_p(k+1)-C_p{\hat{x}}_p(k+1|k)]\end{array}$

其中，和分别表示俯仰通道状态量的预测值和估计值，K_p(k+1)表示俯仰通道的滤波增益阵；

根据俯仰通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_p\left(k\right)\right|\left|\right)$

其中，N表示故障检测窗口的长度，r_p(k)表示俯仰通道第1层状态估计器的测量残差，其形式如下所示：

$r_p\left(k\right)=y_p\left(k\right)-C_p{\hat{x}}_p\left(k\right|k-1)$

基于滚动-偏航通道故障诊断模型的第1层状态估计器为：

${\hat{x}}_r(k+1|k)=A_r{\hat{x}}_r\left(k\right|k)+B_r{u}_r\left(k\right)$

$\begin{array}{c}{\hat{x}}_r(k+1|k+1)=\\ {\hat{x}}_r(k+1|k)+K_r(k+1)[y_r(k+1)-C_r{\hat{x}}_r(k+1|k)]\end{array}$

其中，和分别表示滚动-偏航通道状态量的预测值和估计值，K_r(k+1)表示滚动-偏航通道的滤波增益阵；

根据滚动-偏航通道第1层状态估计器的状态估计结果产生故障检测残差，故障检测残差的计算方式如下所示：

${\mathrm{LD}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}\left(\frac{1}{2}\right|\left|r_r\left(k\right)\right|\left|\right)$

其中，r_p(k)表示滚动-偏航通道第1层状态估计器的测量残差，其形式如下所示：

$r_r\left(k\right)=y_r\left(k\right)-C_r{\hat{x}}_r\left(k\right|k-1)$

(3)如果LD_p(k)≥T_D或LD_r(k)≥T_D，则卫星姿态控制系统中的光学敏感器或惯性敏感器发生了故障，进入步骤(4)；反之则光学敏感器和惯性敏感器均未发生故障；其中，T_D为预设的故障检测阈值；

(4)设计第2层状态估计器；

基于俯仰通道故障诊断模型的第2层状态估计器形式为：

z_p(k+1)＝F_pz_p(k)+G_py_p(k)+T_pB_pu_p(k)

e_p(k)＝K_pz_p(k)+P_py_p(k)

其中，z_p(k)表示俯仰通道第2层状态估计器的状态量，e_p(k)表示俯仰通道第2层状态估计器的输出，F_p、G_p、T_p、K_p和P_p表示俯仰通道第2层状态估计器设计参数矩阵，其中，

F_p＝ρ_pI

I表示单位阵，ρ_p为可调参数且|ρ_p|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_p{A}_p-F_p{T}_p=G_p{C}_p\\ K_p{T}_p+P_p{C}_p=0\\ K_p{T}_p{B}_p=0\end{array}\right.$

根据俯仰通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_p\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_p\left(k\right)$

基于滚动-偏航通道故障诊断模型设计的第2层状态估计器形式为：

z_r(k+1)＝F_rz_r(k)+G_ry_r(k)+T_rB_ru_r(k)

e_r(k)＝K_rz_r(k)+P_ry_r(k)

其中，z_r(k)表示滚动-偏航通道第2层状态估计器的状态量，e_r(k)表示滚动-偏航通道第2层状态估计器的输出，F_r、G_r、T_r、K_r和P_r表示滚动-偏航通道第2层状态估计器设计参数矩阵，其中，

F_r＝ρ_rI

ρ_r为可调参数且|ρ_r|＜1，其它设计参数矩阵通过求解如下所示的矩阵方程得到

$\left\{\begin{array}{c}{T}_r{A}_r-F_r{T}_r=G_r{C}_r\\ K_r{T}_r+P_r{C}_r=0\\ K_r{T}_r{B}_r=0\end{array}\right.$

根据滚动-偏航通道第2层状态估计器的输出产生故障隔离残差，故障隔离残差的计算方式如下所示：

${\mathrm{LI}}_r\left(k\right)=\frac{1}{N}\underset{i=k}{\overset{k+N-1}{\mathrm{\Σ}}}{e}_r\left(k\right)$

(5)如果LI_p(k)≥T₁或LI_r(k)≥T₁，则光学敏感器发生故障，反之，则惯性敏感器发生故障，其中，T₁表示预设的故障隔离阈值；

所述俯仰通道状态量x_P(k)的形式为：

$x_p\left(k\right)=\left[\begin{array}{c}\mathrm{\θ}\left(k\right)\\ d_p\left(k\right)\\ b_p\left(k\right)\\ b_{\mathrm{\θ}}\left(k\right)\end{array}\right]$

其中，θ(k)表示卫星俯仰角，d_p(k)表示俯仰轴陀螺指数漂移误差，b_p(k)表示俯仰轴陀螺常值漂移误差，b_θ(k)表示地球敏感器俯仰角观测量常值偏差；

所述俯仰通道状态矩阵A_p的形式为：

$A_p=\left[\begin{array}{cccc}1& -1& -1& 0\\ 0& 1-{\mathrm{\α}}_p& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right]$

其中，α_p表示俯仰轴陀螺指数偏移误差系数；

所述俯仰通道输入矩阵B_p的形式为：

$B_p=\left[\begin{array}{c}1\\ 0\\ 0\\ 0\end{array}\right].$

所述俯仰通道输入量u_p(k)的形式为：

u_p(k)＝ω₀(k)+g_p(k)

其中，ω₀(k)表示卫星轨道角速率，g_p(k)表示俯仰轴陀螺测量输出；

所述俯仰通道观测量形式为：

$y_p\left(k\right)=\left[\begin{array}{c}{\mathrm{\θ}}_h\left(k\right)\\ {\mathrm{\θ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，θ_h(k)和θ_sh(k)分别表示俯仰轴地球敏感器和太阳敏感器测量输出；

所述俯仰通道观测矩阵形式为： $C_p=\left[\begin{array}{cccc}1& 0& 0& 1\\ 1& 0& 0& 0\end{array}\right];$

所述滚动-偏航通道状态量的形式为：

$x_r\left(k\right)=\left[\begin{array}{c}\mathrm{\φ}\left(k\right)\\ \mathrm{\ψ}\left(k\right)\\ d_r\left(k\right)\\ d_y\left(k\right)\\ b_r\left(k\right)\\ b_y\left(k\right)\\ b_{\mathrm{\φ}}\left(k\right)\end{array}\right]$

其中，φ(k)和ψ(k)分别表示卫星滚动角和偏航角，d_r(k)和d_y(k)分别表示滚动轴和偏航轴陀螺指数漂移误差，b_r(k)和b_y(k)分别表示滚动轴和偏航轴陀螺常值漂移误差，b_φ(k)表示地球敏感器滚动角观测量常值偏差；

所述滚动-偏航通道状态矩阵形式为：

$A_r=\left[\begin{array}{ccccccc}1& {\mathrm{\ω}}_0& -1& 0& -1& 0& 0\\ {-\mathrm{\ω}}_0& 1& 0& -1& 0& -1& 0\\ 0& 0& 1-{\mathrm{\α}}_r& 0& 0& 0& 0\\ 0& 0& 0& 1-{\mathrm{\α}}_y& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 1\end{array}\right]$

其中，α_r和α_y分别表示滚动轴和偏航轴陀螺指数偏移误差系数；

所述滚动-偏航通道输入矩阵形式为：

$B_r=\left[\begin{array}{cc}1& 0\\ 0& 1\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\end{array}\right]$

所述滚动-偏航通道输入量形式为：

$u_r\left(k\right)=\left[\begin{array}{c}{g}_r\left(k\right)\\ g_r\left(k\right)\end{array}\right]$

其中，g_r(k)和g_y(k)分别表示滚动轴和偏航轴陀螺测量输出；

所述滚动-偏航通道观测量形式为：

$y_r\left(k\right)=\left[\begin{array}{c}{\mathrm{\φ}}_h\left(k\right)\\ {\mathrm{\φ}}_{\mathrm{sh}}\left(k\right)\\ {\mathrm{\ψ}}_{\mathrm{sh}}\left(k\right)\end{array}\right]$

其中，φ_h(k)和φ_sh(k)分别表示滚动轴地球敏感器和太阳敏感器测量输出，ψ_sh(k)表示太阳敏感器滚动轴测量输出；

所述滚动-偏航通道观测矩阵形式为：

$C_r=\left[\begin{array}{ccccccc}1& 0& 0& 0& 0& 0& 1\\ 1& 0& 0& 0& 0& 0& 0\\ 0& 1& 0& 0& 0& 0& 0\end{array}\right].$