[发明专利]一种基于捷联惯导与飞行控制系统的大气参数解算方法

摘要

本发明公开了一种基于捷联惯导与飞行控制系统的大气参数解算方法，基于攻角、侧滑角的定义，建立真空速投影到机体系的映射关系，结合高空大气流动速度稳定的特点，利用惯导信息得到大气参数的一步预测模型，根据飞行器气动模型，建立气动参数与大气参数的函数关系，构建起惯导量测值与大气参数的数学关系，通过扩展卡尔曼滤波的方法对大气参数进行求解。该方法为大气参数的测量提供了冗余手段，克服了高超声速、大机动情况下大气参数的测量难点，能够在飞行包线内实时提供大气参数。

主权项

一种基于捷联惯导与飞行控制系统的大气参数解算方法，其特征在于：包括以下步骤：步骤1、初始信息设置：状态量为攻角α、侧滑角β和真空速Vt，设置攻角α的初值为载体俯仰角θ，侧滑角β的初值为零，攻角、侧滑角单位均为rad，真空速Vt的初值为载体速度Vb，单位为ft/s；设置系统噪声方差阵Q、测量噪声方差阵R和初始估计误差协方差阵P；步骤2、状态量变化率解算：根据真空速投影到机体系的映射关系，计算状态量在当前k时刻的变化率，具体步骤包括：步骤201、载入k时刻捷联惯导系统输出的机体系加速度，包括机头方向加速度机身右翼加速度和机身垂向加速度加速度单位均为ft/s2；步骤202、基于真空速投影到机体系的映射关系，进行微分，计算状态量的变化率：α·=cosαkW·k-sinαkU·kcosβkVt,k]]>β·=cosβkV·k-sinβkcosαkU·k-sinβksinαkW·kVt,k]]>V·t=cosαkcosβkU·k+sinβkV·k+cosβksinαkW·k]]>式中，和为k时刻状态量的变化率；αk、βk和Vt,k为k时刻状态量的估计值，初始时刻状态量的估计值取初值；步骤3、状态量一步预测：根据状态量变化率，结合采样步长T，对状态量进行一步预测：α~k+1=αk+α·T]]>β~k+1=βk+β·T]]>V~t,k+1=Vt,k+V·tT]]>式中，和为k时刻状态量的一步预测；步骤4、一步预测均方误差解算：计算状态一步转移矩阵F并离散化，结合k时刻的估计误差协方差阵Pk，解算一步预测均方误差Pk+1,k，具体步骤包括：步骤401、根据k时刻状态量的估计值和机体系加速度，计算状态一步转移矩阵F，F中的元素为：F1,1=-sinαkW·k+cosαkU·kVt,kcosβk]]>F1,2=cosαkW·k-sinαkU·kVt,kcos2βk]]>F1,3=sinαkU·k-cosαkW·kVt,k2cosβk]]>F2,1=sinαksinβkU·k-sinβkcosαkW·kVt,k]]>F2,2=-cosαkcosβkU·k+cosβksinαkW·k+sinβkV·kVt,k]]>F2,3=cosαksinβkU·k+sinβksinαkW·k-cosβkV·kVt,k2]]>F3,1=cosαkcosβkW·k-sinαkcosβkU·k]]>F3,2=cosβkV·k-cosαksinβkU·k-sinβksinαkW·k]]> F3,3＝0步骤402、对F离散化，得到Φk+1,k＝I+FT，I为单位矩阵；步骤403、结合k时刻的估计误差协方差阵Pk，解算一步预测均方误差Pk+1,k：Pk+1,k=Φk+1,kPkΦk+1,kT+Q]]>式中，Pk为k时刻的估计误差协方差阵，初始时刻的估计误差协方差阵为P；步骤5、滤波增益解算：根据捷联惯导系统和飞行控制系统提供的信息，建立气动参数与大气参数的函数关系，解算量测矩阵，具体步骤包括：步骤501、载入下一时刻k+1时刻的捷联惯导系统输出的机体系角速度、飞行控制系统的飞行控制量和空气密度ρ，空气密度的单位为slug/ft3；步骤502、根据气动模型，结合捷联惯导系统和飞行控制系统提供的信息，将气动系数改写为大气参数的函数，解算量测矩阵H，H中元素为：H1,1=d1β~k+1V~t,k+12+d3V~t,k+12+d4V~t,k+1]]>H1,2=d1α~k+1V~t,k+12+d2V~t,k+12]]>H1,3=2d1α~k+1β~k+1V~t,k+1+2d2β~k+1V~t,k+1+2d3α~k+1V~t,k+1+d4α~k+1+d5+2d6V~t,k+1]]>H2,1=d8β~k+1V~t,k+12+d10V~t,k+12+d11V~t,k+1]]>H2,2=d8α~k+1V~t,k+12+d9V~t,k+12]]>H2,3=2d8α~k+1β~k+1V~t,k+1+2d9β~k+1V~t,k+1+2d10α~k+1V~t,k+1+d11α~k+1+d12+2d13V~t,k+1]]>H3,1=d15β~k+1V~t,k+12+d17V~t,k+12+d18V~t,k+1]]>H3,2=d15α~k+1V~t,k+12+d16V~t,k+12]]>H3,3=2d15α~k+1β~k+1V~t,k+1+2d16α~k+1V~t,k+1+2d17α~k+1V~t,k+1+d18α~k+1+d19+2d20V~t,k+1]]>H4,1=d22β~k+1V~t,k+12+d24V~t,k+12+d25V~t,k+1]]>H4,2=d22α~k+1V~t,k+12+d23V~t,k+12]]>H4,3=2d22α~k+1β~k+1V~t,k+1+2d23β~k+1V~t,k+1+2d24α~k+1V~t,k+1+d25α~k+1+d26+2d27V~t,k+1]]>H5,1=2d29α~k+1β~k+1Vt,k+12+2d30α~k+1Vt,k+12+d31β~k+1Vt,k+12+d33Vt,k+12+d34Vt,k+1]]>H5,2=d29α~k+12V~t,k+12+d31α~k+1V~t,k+12+d32V~t,k+12]]>H5,3=2d29α~k+12β~k+1V~t,k+1+2d30α~k+12V~t,k+1+2d31α~k+1β~k+1V~t,k+1+2d32β~k+1V~t,k+1+2d33α~k+1V~t,k+1+d34α~k+1+d35+2d36V~t,k+1]]>H6,1=d38β~k+1V~t,k+12+d40V~t,k+12+d41V~t,k+1]]>H6,2=d38α~k+1V~t,k+12+d39V~t,k+12]]>H6,3=2d38α~k+1β~k+1V~t,k+1+2d39β~k+1V~t,k+1+2d40α~k+1V~t,k+1+d41α~k+1+d42+2d43V~t,k+1]]> 式中，di为函数中的常数系数，由气动系数的函数系数、载体的速度、角速度、姿态以及转动惯量解算获得，i＝1～43；步骤503、解算滤波增益步骤6、估计均方误差解算：根据滤波增益Kk+1、量测矩阵H、一步预测均方误差Pk+1,k和测量噪声方差阵R，解算估计均方误差Pk+1：Pk+1=(I-Kk+1H)Pk+1,k(I-Kk+1H)T+Kk+1RKk+1T;]]>步骤7、状态估计：结合状态量的一步预测，根据载体的气动模型，解算k+1时刻的观测量估计值，从而实现状态量的估计，即大气参数的估计，具体步骤包括：步骤701、载入k+1时刻捷联惯导系统输出的机体系角加速度，包括机头方向角加速度机身右翼角加速度和机身垂向角加速度角加速度单位均为rad/s2；载入k+1时刻的飞行控制量、机体系加速度步骤702、将状态量在k时刻的一步预测、k+1时刻的飞行控制量、k+1时刻的机体系加速度和k+1时刻的机体系角加速度代入到气动模型中，解算观测量估计值：U·~k+1=d1α~k+1β~k+1V~t,k+12+d2β~k+1V~t,k+12+d3α~k+1V~t,k+12+d4α~k+1V~t,k+1+d5V~t,k+1+d6V~t,k+12+d7]]>V·~k+1=d8α~k+1β~k+1V~t,k+12+d9β~k+1V~t,k+12+d10α~k+1V~t,k+12+d11α~k+1V~t,k+1+d12V~t,k+1+d13V~t,k+12+d14]]>W·~k+1=d15α~k+1β~k+1V~t,k+12+d16β~k+1V~t,k+12+d17α~k+1V~t,k+12+d18α~k+1V~t,k+1+d19V~t,k+1+d20V~t,k+12+d21]]>P·~k+1=d22α~k+1β~k+1V~t,k+12+d23β~k+1V~t,k+12+d24α~k+1V~t,k+12+d25α~k+1V~t,k+1+d26V~t,k+1+d27V~t,k+12+d28]]>Q·~k+1=d29α~k+12β~k+1V~t,k+12+d30α~k+12V~t,k+12+d31α~k+1β~k+1V~t,k+12+d32β~k+1Vt,k+12+d33α~k+1V~t,k+12+d34α~k+1V~t,k+1+d35V~t,k+1+d36V~t,k+12+d37]]>R·~k+1=d38α~k+1β~k+1V~t,k+12+d39β~k+1V~t,k+12+d40α~k+1V~t,k+12+d41α~k+1V~t,k+1+d42V~t,k+1+d43V~t,k+12+d44]]>式中，d44为函数中的常数，由气动系数的函数系数、载体的速度、角速度、姿态以及转动惯量解算获得；步骤703、记k+1时刻的观测量为k+1时刻的观测量估计值为状态量在k时刻的一步预测为状态量在k+1时刻的一步预测为则有状态估计值在获得状态估计值后，返回步骤2，继续进行后续的大气参数解算。