[发明专利]一种挠性卫星模态参数在轨辨识方法

摘要

一种挠性卫星模态参数在轨辨识方法，本发明涉及卫星模态参数在轨辨识领域，本发明要解决无法建立精确的动力学模型，频域法时效性低，难以辨识密频模态参数及获取阻尼信息使得估算、综合的参数存在误差，辨识精度降低以及时域法无法准确对模型定阶的问题而提出的一种挠性卫星模态参数在轨辨识方法。该方法是通过1、采集力矩和角速度信息；2、确定离散状态空间系统矩阵A、B、C和D；3、确定辨识矩阵A、B、C和D；4、建立动力学和帆板振荡方程；5、获得模态参数与力矩到角速度的传递函数；6、将矩阵A、B、C和D转换为传递函数等步骤实现的。本发明应用于卫星模态参数在轨辨识领域。

主权项

一种挠性卫星模态参数在轨辨识方法，其特征在于：一种挠性卫星模态参数在轨辨识方法具体是按照以下步骤进行的：步骤一、在外界诸多干扰不影响挠性卫星姿态稳定的前提下，于每个采样周期T平均准确的采集到执行机构施加到挠性卫星体上的力矩和挠性卫星体相对惯性坐标系的角速度信息；步骤二、将采集到的力矩和角速度信息代入子空间辨识算法，确定离散状态空间系统矩阵A和C并且计算出离散状态空间系统矩阵B和D；步骤三、假设当前所需辨识系统的系统状态方程，确定离散状态空间系统矩阵A、B、C和D；步骤四、基于挠性卫星姿态角的小角度假设，建立挠性卫星动力学方程和帆板振荡方程；步骤五、将建立的挠性卫星动力学方程和帆板振荡方程解耦，通过拉式变换，获得模态参数与力矩到角速度的传递函数；步骤六、将步骤三得到的离散状态空间系统矩阵A、B、C和D转换为步骤五获得的模态参数与力矩到角速度的传递函数；即完成了一种挠性卫星模态参数在轨辨识方法；步骤二中将采集到的力矩和角速度信息代入子空间辨识算法，确定离散状态空间系统矩阵A和C并且计算出离散状态空间系统矩阵B和D的实现步骤为：(一)通过输入输出数据构造Hankel块矩阵；其中输入块Hankel矩阵由输入数据构造，具体构成形式：其中，(1)i表示行块数，每一个行块含有m行，U0|2i‑1有2mi行；j为d‑2i+1；(2)输入块Hankel矩阵U0|2i‑1的脚标由该输入块Hankel的第一列中第一个行块u0和最后一个行块u2i‑1确定，输入块式中“p”表示“过去”，“f”表示“将来”；那么矩阵Up表示由过去的输入数据构成的过去输入块Hankel矩阵，矩阵Uf表示由将来的输入数据构成的将来输入块Hankel矩阵；将矩阵U0|2i‑1分为两个矩阵块U0|i‑1,Ui|2i‑1，并将U0|i‑1定义为矩阵Up，将Ui|2i‑1定义为矩阵Uf；输出块Hankel矩阵Y0|2i‑1，Yp，Yf具体定义形式如下：其中，矩阵Yp表示由过去的输出数据构成的输出块Hankel矩阵，矩阵Yf表示由将来的输出数据构成的输出块Hankel矩阵，输出块Hankel矩阵Y0|2i‑1由该输出块Hankel的第一列中第一个行块y0和最后一个行块y2i‑1确定；包含过去输入Up和过去输出Yp的过去输入输出Hankel块矩阵Wp定义：Wp=UPYP=U0|i-1Y0|i-1---(2)]]>状态序列Xi定义为：其中，xi表示状态序列的第一个矩阵块,xi+j‑1表示状态序列的最后一个矩阵块，状态序列为n行j列，表示状态序列，状态序列Xi表示维数为n行j列实数集矩阵；其中i表示状态序列的第一个元素的下标为每个矩阵行块数；同时有：Xp＝X0,Xf＝Xi  (4)其中，Xp表示过去的实数集构成的矩阵，X0表示过去的实数集构成的矩阵第一个元素构成的矩阵，Xf表示将来的实数集构成的矩阵；计算将来输出块Hankel矩阵Yf的行空间沿着将来输入块Hankel矩阵Uf包含过去输入输出的Hankel矩阵行空间上的斜投影Oi；Oi=Yf/UfWp---(5)]]>公式中：Uf：定义为将来输入块Hankel矩阵；Yf：定义为将来输出块Hankel矩阵；Wp：包含有过去输入输出的Hankel矩阵；其中，正交投影为：∏B表示将一个矩阵的行空间投影到离散状态空间系统矩阵的行空间的几何操作：其中：表示矩阵的伪逆，矩阵表示矩阵BBT的伪逆矩阵；表示将一个矩阵的行空间投影到矩阵的行空间的正交空间的几何操作，并且∏B和存在式(7)的关系：ΠB⊥=Ij-ΠB---(7)]]>其中：Ij：为j阶单位方阵A/B表示将离散状态空间系统矩阵的行空间投影到B上的几何操作：那么A/B⊥=AΠB⊥---(9)]]>矩阵的行空间沿着的行空间在的行空间上的斜投影定义为：A/BC=ACTBTCCTCBTBCTBBTfirstrcolumnsC---(10)]]>其中，first r columns：表示矩阵的前r列；将A的行空间向B和C的联合行空间正交投影，即：A/BC=A/CB---(11)]]>(二)计算加权斜投影SVD：对式(5)进行奇异值分解，分解为：W1OiW2=USVT=U1U2S1000V1TV2T=U1S1V1T---(12)]]>W1∈Rli×li：自定义权值矩阵，为满秩矩阵；W2∈Rj×j：自定义权值矩阵，S为奇异值分解过程中的矩阵，U,V均为酉矩阵，酉矩阵U的列向量被称为Oi的左奇异向量，将U从前r列处分块为U＝(U1 U2)，酉矩阵V的列向量称为Oi的右奇异向量，V的前r列是OiHOi的r个非零特征值所对应的特征向量，则V＝(V1,V2)，且满足：rank(Wp)＝rank(WP W2)；则矩阵Oi满足：Oi＝ΓiXf  (13)其中，Γi：扩展可观测矩阵，满秩；表示为A,C：均为离散状态空间系统矩阵，而且离散状态空间系统矩阵A,C均可观，i：表示每个矩阵行块数，n为系统阶数，l表示输出量的个数，Xf为状态序列，辨识的线性系统阶数n等于式(12)中非零奇异值的个数；式(12)中奇异值分解过程中的矩阵S的奇异值确定系统阶数n，并获得扩展可观矩阵为：状态序列Xf中由W2的列空间决定的部分表示为：XfW2=A-1S11/2V1T---(15)]]>其中，A为离散状态空间系统矩阵；状态序列Xf等于：(三)由矩阵U,S,W确定扩展可观测矩阵Γi和扩展可观测矩阵Γi行空间的正交补矩阵Γi=W1-1U1S11/2---(17)]]>Γi⊥=U2TW1---(18)]]>其中，Γi：为扩展观测矩阵；扩展可观测矩阵Γi行空间的正交补矩阵；(四)由利用输入输出矩阵获得相关扩展可观测矩阵Γi和扩展可观测矩阵Γi行空间的正交补矩阵确定离散状态空间系统矩阵A，C；离散状态空间系统矩阵A由扩展可观测矩阵Γi的移位后矩阵计算得到：Γi‾A=Γi‾---(19)]]>那么得到：C=Γifirstlrow---(21)]]>其中，表示扩展可观测矩阵Γi去掉前l行的矩阵；Γi：表示扩展可观测矩阵Γi去掉后l行的矩阵；表示矩阵的伪逆，即表示矩阵Γi的伪逆；Γi的前l行；(五)辨识系统矩阵矩阵输入输出方程定义为：Yp＝ΓiXp+HiUp  (22)Yf＝ΓiXf+HiUf  (23)Xf＝AiXp+ΔiUp  (24)其中，Γi：表示系统扩展可观测矩阵，定义为：Hi：表示三角形式的Toeplitz矩阵，定义为：Δi：表示系统扩展控矩阵，定义为：其中，n表示确定系统的阶数；m表示输入的维数；公式(23)左侧乘以得到公式(25)：Γi⊥Yf=Γi⊥Hi.Uf---(25)]]>在式(25)的右侧乘以得到公式(26)：将矩阵M代替式(26)中即：用矩阵L代替式(26)右侧的扩展可观测矩阵即L=L1L2...Li=Γi⊥---(28)]]>其中，为将来输入块Hankel矩阵Uf的伪逆；那么式(26)表示为：M1M2...Mi=L1L2...Li×D00...0CBD0...0CABCBD...0...............CAi-2BCAi-3BCAi-4B...D---(29)]]>写成如下形式：离散状态空间系统矩阵B和D由式(29)(30)计算得到。