CN115615261B - 大长细比火箭弹弹性辨识与视线角速率提取融合方法 - Google Patents

Abstract

本发明公开了一种大长细比火箭弹弹性辨识与视线角速率提取融合方法，包括以下步骤：S1、以视线角、视线角速率、姿态角、转动角速度、火箭弹振动方程广义坐标及其广义坐标的变化率作为状态变量，建立大长细比火箭弹视线角速率状态模型；S2、对火箭弹体视线高低角、火箭弹体视线方位角进行修正；S3、以火箭弹体视线高低角、火箭弹体视线方位角以及火箭弹弹体三个轴向的角速度作为观测变量，建立大长细比火箭弹视线角速率观测模型；S4、基于视线角速率状态模型和观测模型，通过滤波估计获得视线角速率。本发明提供的方法，实现了弹性形变情况下视线角速率的高精度提取，且提取准确度不会受到信号传递时延影响。

Description

大长细比火箭弹弹性辨识与视线角速率提取融合方法

技术领域

本发明涉及一种大长细比火箭弹弹性辨识与视线角速率提取融合方法，属于飞行器控制领域。

背景技术

大长细比构型并配备捷联导引头的制导火箭弹能够有效提高火箭弹的射程 和精度，但是该构型的火箭弹在飞行过程中受到气动弹性力等因素的影响会产生 周期性的弹性形变，尤其是在被动飞行阶段，随着固体燃料的消耗，使得火箭弹 弹体刚度变弱，发生弹性形变的幅度大大增加。对于火箭弹上搭载的捷联导引头 来说，由于敏感元件与弹体固连，无法隔离弹性形变，导致弹性形变混入导引头的目标观测信息与陀螺仪的姿态信息。

传统视线角速率获取方法是基于刚体假设提出的，用于为导弹提供准确的视 线角速率及视线角，但是对于配备捷联导引头的大长细比火箭弹，用来提供姿态 解耦信息的角速率陀螺以及用来提供目标角度信息的导引头将受到火箭弹弹性 形变的干扰，导致传统方法提取到的视线角速率具有较大偏差。

实际上，大长细比火箭弹在整个无控飞行过程中，纵向弹性形变量则始终存 在，并逐渐增大，如图2所示，尤其是在弹道末段由于飞行速度快、弹道高度较 低、动压变大，导致弹性形变幅度不断增大的变化趋势，受到火箭弹的弹性形变 的影响，弹体俯仰方向的速度、角速度以及攻角在都产生了一定幅度的振动，进 而影响到火箭弹姿态传感器的输出。如果不进行处理，会造成火箭弹制导和姿态 控制等环节产生较大误差。

针对这种情况，目前主流研究方向致力于对弹性形变进行主动抑制，来减弱 甚至消除弹性形变对弹体传感器的影响，例如文献PAN C L,RONG J L,XU T F, et al.Novelapproach for active vibration control of a flexible missile[J].DefenceTechnology,2020,16(4):956-967.，通过给弹体加入一个环形弹性抑制机构，利用 遗传算法设计PID控制器实现对弹性形变的主动抑制，然而该方法最多只能实 现约40％的抑制效果，又例如文献张博伦,周荻.引入角加速度测量的柔性飞行 器姿态控制方法[J].兵工学报,2020,41(11):2225-2233.通过控制舵机产生与弹 性形变相反方向的控制力来维持弹体的稳定，然而]通过数值仿真发现，当弹性 形变引起的扰动混入加速度和角速度传感器，会引起闭环反馈系统严重发散，并 且大长细比火箭弹的结构频率与刚体运动频率更加接近，很难从刚体运动中分离 弹性形变。

由于上述原因，有必要提出一种能够解决上述问题的大长细比火箭弹弹性辨识与视线角速率提取融合方法。

发明内容

为了克服上述问题，本发明人进行了深入研究，设计出一种大长细比火箭弹弹性辨识与视线角速率提取融合方法，包括以下步骤：

S1、以视线角、视线角速率、姿态角、角速度、火箭弹振动方程广义坐标及其广 义坐标的变化率作为状态变量，建立大长细比火箭弹视线角速率状态模型；

S2、对火箭弹体视线高低角、火箭弹体视线方位角进行修正；

S3、以火箭弹体视线高低角、火箭弹体视线方位角以及火箭弹弹体三个轴向的角速度作为观测变量，建立大长细比火箭弹视线角速率观测模型；

S4、基于视线角速率状态模型和观测模型，通过滤波估计获得视线角速率。

进一步地，在S1中，所述大长细比火箭弹视线角速率状态模型表示为：

其中，状态变量[x₁,x₂,x₃,x₄,x₅,x₆,x₇,x₈,x₉,x₁₀,x₁₁,x₁₂,x₁₃,x₁₄]^T设置为：

q_z表示火箭弹视线高低角，q_y表示火箭弹视线方位角；θ、ψ、γ为地面固连坐标系下 火箭弹的姿态角；ω_x、ω_y、ω_z为瞬态坐标系下火箭弹的角速度；η_y、η_z为火箭弹振 动方程拉格朗日函数第1阶模态下的广义坐标；r表示火箭弹与目标之间的相对距离；/>

表示火箭弹在体视线坐标系内垂直方向加速度分量，/>

表示火箭弹在体视线坐标系 内水平加速度分量；M_x、M_y、M_z为火箭弹质心处气动弹性力矩；J_x为火箭弹极转动 惯量，J_y与J_z为赤道转动惯量；N_y、N_z分别为瞬态坐标系下垂直方向和侧向广义力； M₁为瞬态坐标系下火箭弹动能的拉格朗日函数第1阶模态的质量；μ₁为瞬态坐标系下 火箭弹势能的拉格朗日函数第1阶模态的临界阻尼系数，ω₁为瞬态坐标系下火箭弹势 能的拉格朗日函数第1阶模态的固有频率。

进一步地，设置火箭弹相对地面坐标系姿态变化的运动学模型对参数θ、ψ、γ、ω_x、ω_y、ω_z进行约束，所述火箭弹相对地面坐标系姿态变化的运动学模型表示为：

其中，t表示时间。

进一步地，设置绕质心动力学模型对参数ω_x、ω_y、ω_z进行约束，所述绕质心动 力学模型表示为：

其中，

b₁＝-σ₁η_z,1；

c₁＝-σ₁η_y,1；

a₂＝-σ₁η_y,1；

a₃＝-σ₁η_z,1

η_y,1表示垂直方向弹性形变的广义坐标，η_z,1表示侧向弹性形变的广义坐标，m为火箭弹质量，σ₁为常数，表征弹性形变与刚体旋转耦合引起的动能增量，λ₁为常数， 表征引起火箭弹质心位移大小。

进一步地，设置瞬态坐标系下的大长细比火箭弹的动力学模型对参数M_x、M_y、 M_z、η_y、η_z、N_y、N_z、ω_x、ω_y、ω_z进行约束，所述瞬态坐标系下的大长细比火箭弹 的动力学模型表示为：

其中，V_x、V_y和V_z分别为火箭弹在地面固连坐标系下的速度分量，L为拉格朗日 函数，表示为L＝T-P，T为火箭弹动能，P为火箭弹势能，D为阻尼耗散能，F_x、F_y、 F_z为火箭弹所受合力在瞬态坐标系下的分量。

进一步地，所述拉格朗日函数第1阶模态下的振动方程表示为：

其中，λ₁为引起火箭弹质心位移大小的常数，σ₁为弹性形变与刚体旋转耦合引起的动能增量有关的常数，M₁为为瞬态坐标系下火箭弹动能的拉格朗日函数第1阶模态 质量，通过火箭弹的动能T获得。

进一步地，所述火箭弹的动能T表示为：

其中，

M_i＝∫_mΦ_iΦ_idm，λ_i＝∫_mΦ_idm，

当n＝1时，通过火箭弹的动能T即可获得λ₁、σ₁、M₁。

在一个优选的实施方式中，在S2中，以导引头位置处的弹性形变角对导引头体视线角进行修正，表示为：

q_αm＝q_α-Φ′₁(x_head)η_y，q_βm＝q_β-Φ′₁(x_head)η_z，

q_αm、q_βm为修正后的导引头体视线角，q_α、q_β为刚性运动的体视线角，由导引头 输出，下标α表示体视线坐标系下导引头垂直方向分量，下标β表示体视线坐标系下 导引头水平方向分量，Φ′₁(x_head)η_y、Φ′₁(x_head)η_z为导引头位置处的弹性形变角，x_head表示 导引头位置，

设置全弹气动力和力矩模型对导引头位置处的弹性形变角Φ′₁(x_head)η_y、Φ′₁(x_head)η_z进行约束；

所述全弹气动力和力矩模型表示为：

其中，

为火箭弹轴段微元x_b的升力系数导数、/>

为火箭弹轴段微元x_b的侧向力系数导数，Q为动压，S为火箭弹截面积，l表示火箭弹弹长，x_G为火箭弹末端 到火箭弹质心的轴向距离；α(x_b,t)为在弹性形变影响下，t时刻在弹体轴向任一微元x_b的局部攻角，β(x_b,t)为在弹性形变影响下，t时刻在弹体轴向任一微元x_b的侧滑角， Φ′₁(x)表示Φ₁(x)对x的导数。

在一个优选的实施方式中，在S3中，所述大长细比火箭弹视线角速率观测模型表示为：

其中，观测变量[z₁,z₂,z₃,z₄z₅]^T设置为：

[z₁,z₂,z₃,z₄z₅]^T＝[q_αm,q_βm,ω_xm,ω_ym,ω_zm]^T；

x_imu表示火箭弹陀螺仪安装位置，ω_xm,ω_ym,ω_zm为火箭弹耦合弹性形变后的陀螺仪角 速度的输出，a、b、c为中间量。

在一个优选的实施方式中，S4中，基于视线角速率状态模型和观测模型，通过无迹卡尔曼滤波法进行滤波估计，获得视线角速率。

本发明所具有的有益效果包括：

(1)实现了弹性形变情况下视线角速率的高精度提取；

(2)经过仿真研究表明，提取得到的视线角速率误差最大值在0.2°以内；

(3)将弹性辨识和视线角速率提取融合，不会受到信号传递时延影响。

附图说明

图1出根据本发明一种优选实施方式的无位置传感器永磁同步电机全速域带速重投的方法整体流程示意图；图2示出根据本发明一种优选实施方式的无位置传感器永 磁同步电机全速域带速重投的方法中单次短路电机三相绕组电流矢量示意图；图3示 出了实施例1与对比例1中视线高低角估计曲线；图4示出了实施例1与对比例1中视线 方位角估计曲线；图5示出了实施例1与对比例1中视线高低角估计误差曲线；图6示出了实施例1与对比例1中视线方位角估计误差曲线；图7示出了实施例1与对比例2中视线 高低角估计曲线；图8示出了实施例1与对比例2中视线方位角估计曲线；图9示出了实施 例1与对比例2中视线高低角估计误差曲线。

具体实施方式

下面通过附图和实施例对本发明进一步详细说明。通过这些说明，本发明的特点和优点将变得更为清楚明确。

在这里专用的词“示例性”意为“用作例子、实施例或说明性”。这里作为“示例 性”所说明的任何实施例不必解释为优于或好于其它实施例。尽管在附图中示出了实 施例的各种方面，但是除非特别指出，不必按比例绘制附图。

由于大长细比火箭弹在飞行末段高度、速度变化较大，导致火箭弹所受的气动力急剧变化，会使弹体受力产生弹性形变，这些形变会混杂在导引头、陀螺仪角速度输 出信号中，如果不加处理，将对火箭弹的制导系统和控制系统产生不利影响。特别是 对于视线角速率的提取，弹体的弹性形变会影响到导引头体视线角输出以及导引头姿态解耦的计算，因此，需要对大长细比火箭弹视线角速率进行有效的获取。

本发明提供的一种大长细比火箭弹弹性辨识与视线角速率提取融合方法，所述大长细比火箭是指弹体长度与弹体直径之比大于24的火箭弹。

进一步地，所述大长细比火箭弹弹性辨识与视线角速率提取融合方法包括以下步骤，如图1所示：

S1、以视线角、视线角速率、姿态角、转动角速度、火箭弹振动方程广义坐标及 其广义坐标的变化率作为状态变量，建立大长细比火箭弹视线角速率状态模型；

S2、对火箭弹体视线高低角、火箭弹体视线方位角进行修正；

S3、以火箭弹体视线高低角、火箭弹体视线方位角以及火箭弹弹体三个轴向的角速度作为观测变量，建立大长细比火箭弹视线角速率观测模型；

S4、基于视线角速率状态模型和观测模型，通过滤波估计获得视线角速率。

由于火箭弹的弹性变形量会影响视线角速率提取的准确率，发明人发现，虽然可以通过弹性辨识过程获取火箭弹的弹性变形量，再利用获取的弹性变形量以及将弹性 变形量进行火箭弹视线角速率提取的过程会出现变形量传递时延，这个时延将导致视 线角速率提取过程中引入的弹性形变不能与实际弹性形变相匹配，导致获得的角速率准确率降低。

在本发明S1中，通过将火箭弹的弹性辨识过程和视线角速率提取过程进行合并研究，使辨识和提取可以在同一个采样周期完成，从而避免采样信号传递出现时延，提 高提取的结果的准确性。

具体地，通过以视线角、视线角速率、姿态角、转动角速度、火箭弹振动方程广 义坐标及其广义坐标的变化率作为状态变量，其中，视线角、视线角速率与实现角速 率获取直接关联，姿态角、转动角速度、火箭弹振动方程广义坐标及其广义坐标的变 化率与火箭弹的弹性变形量直接关联，通过将这些量融合为一体，使得建立的大长细比火箭弹视线角速率状态模型将弹性辨识过程和视线角速率获取过程融合，从而解决 上述问题。

进一步地，在S1中，所述大长细比火箭弹视线角速率状态模型基于拉格朗日函数构建，表示为：

其中，状态变量[x₁,x₂,x₃,x₄,x₅,x₆,x₇,x₈,x₉,x₁₀,x₁₁,x₁₂,x₁₃,x₁₄]^T设置为：

q_z表示火箭弹视线高低角，q_y表示火箭弹视线方位角；θ、ψ、γ为地面固连坐标 系下火箭弹的姿态角；ω_x、ω_y、ω_z为瞬态坐标系下火箭弹的角速度；η_y、η_z为火箭 弹振动方程拉格朗日函数第1阶模态下的广义坐标；r表示火箭弹与目标之间的相对 距离；

表示火箭弹在体视线坐标系内垂直方向加速度分量，/>

表示火箭弹在体视线 坐标系内水平加速度分量；M_x、M_y、M_z为火箭弹质心处气动弹性力矩；J_x为火箭弹 极转动惯量，J_y与J_z为赤道转动惯量；N_y、N_z分别为瞬态坐标系下垂直方向和侧向 广义力，通过全弹气动力和力矩模型获得；M₁为瞬态坐标系下火箭弹动能的拉格朗 日函数第1阶模态的质量；

μ₁为瞬态坐标系下火箭弹势能的拉格朗日函数第1阶模态的临界阻尼系数，一般取0.05～0.1，ω₁为瞬态坐标系下火箭弹势能的拉格朗日函数第1阶模态的固有频率。

所述地面固连坐标系是以发射点为原点A，x轴位于水平面内，指向目标为正； y轴位于铅垂面内，竖直向上为正的坐标系，z轴和x、y轴垂直并构成右手定则的坐标系；

所述瞬态坐标系是指，以弹性飞行器质心为原点O_b，O_bx_b轴与未变形弹体纵轴平行，指向未变形弹体头部为正；O_by_b轴与未变形弹体纵向对称面平行，并垂直于O_bx_b轴，向上为正；O_bz_b轴垂直于x_bO_by_b平面，方向由右手定则确定的坐标系。

在本发明中，所述大长细比火箭弹视线角速率状态模型基于拉格朗日函数构建，利用达朗贝尔原理，即在虚位移上，主动力所做虚功与惯性力所做虚功之和为零，实 现了从系统的能量与功的角度来描述耦合了弹性形变的大长细比火箭弹运动过程。

与传统表征相同，在本发明中，参数上方标记一点表示对参数求一阶导数，标记两点表示二阶导数。

进一步地，设置火箭弹相对地面坐标系姿态变化的运动学模型对参数θ、ψ、γ、ω_x、ω_y、ω_z进行约束，所述火箭弹相对地面坐标系姿态变化的运动学模型表示为：

其中，t表示时间。

进一步地，设置绕质心动力学模型对参数ω_x、ω_y、ω_z进行约束，所述绕质心动 力学模型表示为：

其中，

b₁＝-σ₁η_z,1；

c₁＝-σ₁η_y,1；

a₂＝-σ₁η_y,1；

a₃＝-σ₁η_z,1

η_y,1表示垂直方向弹性形变的广义坐标，η_z,1表示侧向弹性形变的广义坐标，m为火箭弹质量，σ₁为常数，表征弹性形变与刚体旋转耦合引起的动能增量，λ₁为常数， 表征引起火箭弹质心位移大小。

进一步地，设置瞬态坐标系下的大长细比火箭弹的动力学模型对参数M_x、M_y、M_z、η_y、η_z、N_y、N_z、ω_x、ω_y、ω_z进行约束，所述瞬态坐标系下的大长细比火箭弹 的动力学模型表示为：

其中，V_x、V_y和V_z分别为火箭弹在地面固连坐标系下的速度分量，L为拉格朗日 函数，表示为L＝T-P，T为火箭弹动能，P为火箭弹势能，D为阻尼耗散能，F_x、F_y、 F_z为火箭弹所受合力在瞬态坐标系下的分量。

进一步地，火箭弹势能P＝P_e+P_g，其中，P_e为弹性势能，P_g为重力势能，其表达 式分别为：

P_g＝-∫_mR₀gdm＝mgy

其中，R₀为地面固连坐标系到大长细比火箭弹质心的位置矢量；

阻尼耗散能D表示为：

其中，μ_i为第i阶模态的临界阻尼系数，一般取0.05～0.1，ω_i为第i阶模态的固有频率。

进一步地，火箭弹所受合力包括：发动机推力T_b、弹性气动力F_e和重力F_G；优选 地，合力在瞬态坐标系下的分量表示为：

其中，

为发动机喷口处弹性形变量组成的坐标转换矩阵，T_m为火箭发动机喷 口处推力，/>

为地面坐标系到瞬态坐标系的转换矩阵；m为火箭弹质量，g为重力加 速度，

具体地，

根据本发明，所述拉格朗日函数第1阶模态下的振动方程表示为：

其中，λ₁为引起火箭弹质心位移大小的常数，σ₁为弹性形变与刚体旋转耦合引起的动能增量有关的常数，通过火箭弹的动能T获得。

所述火箭弹的动能T可以表示为：

其中，

M_i＝∫_mΦ_iΦ_idm，λ_i＝∫_mΦ_idm，Φ_i为第i阶模态的振型函数，

sinh为双曲正弦函数，cosh为双曲余弦函数，

当n＝1时，通过火箭弹的动能T即可获得λ₁、σ₁、M₁。

在S2中，以导引头位置处的弹性形变角对导引头体视线角进行修正，表示为：

q_αm＝q_α-Φ′₁(x_head)η_y，q_βm＝q_β-Φ′₁(x_head)η_z，

q_αm为修正后的导引头垂直方向体视线角、q_βm为修正后的导引头水平方向体视线角，

q_α为体视线坐标系下刚性运动导引头垂直方向的体视线角、q_β为体视线坐标系下刚性运动导引头水平方向的体视线角，均由导引头输出，

Φ′₁(x_head)η_y为导引头位置处的弹性形变角的垂直分量、Φ′₁(x_head)η_z为导引头位置处的 弹性形变角的侧向分量，x_head表示导引头位置，

设置全弹气动力和力矩模型对导引头位置处的弹性形变角Φ′₁(x_head)η_y、Φ′₁(x_head)η_z进行约束；

所述全弹气动力和力矩模型表示为：

其中，x_b为火箭弹轴段微元，

为火箭弹轴段微元x_b的升力系数导数、/>

为火箭弹轴段微元x_b的侧向力系数导数，Q为动压，S为火箭弹截面积，l表示火箭弹 弹长，x_G为火箭弹末端到火箭弹质心的轴向距离；α(x_b,t)为在弹性形变影响下，t时 刻在弹体轴向任一微元x_b的局部攻角，β(x_b,t)为在弹性形变影响下，t时刻在弹体轴向 任一微元x_b的侧滑角，Φ′₁(x)表示Φ₁(x)对x的导数。

进一步地，α(x_b,t)和β(x_b,t)可以表示为：

α₀、β₀分别为刚体攻角、侧滑角；

分别为弹体弹性变形对时间的导数；d′_y＝Φ′₁(x_b)η_y(t)、d′_z＝Φ′₁(x_b)η_z(t)分别为弹体弹性变形对轴向坐标的 导数。

在S3中，所述大长细比火箭弹视线角速率观测模型表示为：

其中，观测变量[z₁,z₂,z₃,z₄z₅]^T设置为：

[z₁,z₂,z₃,z₄z₅]^T＝[q_αm,q_βm,ω_xm,ω_ym,ω_zm]^T；

x_imu表示火箭弹陀螺仪安装位置，ω_xm,ω_ym,ω_zm为火箭弹耦合弹性形变后的陀螺仪角 速度的输出，a、b、c为中间量。

根据本发明一个优选的实施方式，S4中，基于视线角速率状态模型和观测模型，通过无迹卡尔曼滤波法进行滤波估计，获得视线角速率。

无迹卡尔曼滤波法(Unscented Kalman Filter,UKF)是一种经典的对系统状 态进行最优估计的方法，在本发明中不做赘述。

具体地，对视线角速率提取状态模型和观测模型进行离散化，对离散化后的结果采用无迹卡尔曼滤波(UKF)法进行滤波，无迹卡尔曼滤波法是对非线性函数的概率 密度分布进行近似，而不是对非线性函数进行近似，其选用确定的样本Sigma点来逼 近系统状态的后验概率密度，能较好地描述状态的真实均值和方差，避免了计算Jacobian矩阵，且在高斯条件下可达到三阶Taylor级数展开精度。

无迹卡尔曼滤波法通过在估计点附近进行UT变换，使获得的Sigma点集的状态 均值和协方差与原统计特性匹配，再直接对这些点集进行非线性映射，以近似得到状态概率密度函数，能够以较低的计算复杂度获得较高的估计精度，实现对视线角速率 的高精度提取。

实施例

实施例1

基于装载捷联导引头的某型远程制导火箭弹的无控弹道仿真实验，设置仿真初始参数如下表所示：

通过以下步骤获得火箭弹视线角速率：

S1、以视线角、视线角速率、姿态角、转动角速度、火箭弹振动方程广义坐标及 其广义坐标的变化率作为状态变量，建立大长细比火箭弹视线角速率状态模型；

S2、对火箭弹体视线高低角、火箭弹体视线方位角进行修正；

S3、以火箭弹体视线高低角、火箭弹体视线方位角以及火箭弹弹体三个轴向的角速度作为观测变量，建立大长细比火箭弹视线角速率观测模型；

S4、基于视线角速率状态模型和观测模型，通过滤波估计获得视线角速率。

在S1中，所述大长细比火箭弹视线角速率状态模型表示为：

其中，状态变量[x₁,x₂,x₃,x₄,x₅,x₆,x₇,x₈,x₉,x₁₀,x₁₁,x₁₂,x₁₃,x₁₄]^T设置为： [-0.8464,-0.0025,0,0,-0.7752,0,0,0,0,-0.0072,-0.0008,-0.0133,0,0]^T；μ₁取0.05。

设置火箭弹相对地面坐标系姿态变化的运动学模型对参数θ、ψ、γ、ω_x、ω_y、 ω_z进行约束，所述火箭弹相对地面坐标系姿态变化的运动学模型表示为：

设置绕质心动力学模型对参数ω_x、ω_y、ω_z进行约束，所述绕质心动力学模型表示为：

其中，

b₁＝-σ₁η_z,1；

c₁＝-σ₁η_y,1；

a₂＝-σ₁η_y,1；

a₃＝-σ₁η_z,1

设置瞬态坐标系下的大长细比火箭弹的动力学模型对参数M_x、M_y、M_z、η_y、η_z、 N_y、N_z、ω_x、ω_y、ω_z进行约束，所述瞬态坐标系下的大长细比火箭弹的动力学模型 表示为：

/>

所述拉格朗日函数第1阶模态下的振动方程表示为：

火箭弹的动能T表示为：

其中，

其中截面积S＝0.07068m²，弹长l＝7.6m，

当n＝1时，通过火箭弹的动能T即可获得λ₁、σ₁、M₁。

在S2中，以导引头位置处的弹性形变角对导引头体视线角进行修正，表示为：

q_αm＝q_α-Φ′₁(x_head)η_y，q_βm＝q_β-Φ′₁(x_head)η_z，

设置全弹气动力和力矩模型对导引头位置处的弹性形变角Φ′₁(x_head)η_y、Φ′₁(x_head)η_z进行约束；

所述全弹气动力和力矩模型表示为：

在S3中，所述大长细比火箭弹视线角速率观测模型表示为：

其中，观测变量[z₁,z₂,z₃,z₄z₅]^T设置为：

[z₁,z₂,z₃,z₄z₅]^T＝[q_αm,q_βm,ω_xm,ω_ym,ω_zm]^T；

S4中，基于视线角速率状态模型和观测模型，通过无迹卡尔曼滤波法进行滤波估计，获得火箭弹视线角速率。

对比例1

进行与实施例1相同的实验，区别在于，采用文章ZHANG D,SONG J M,ZHAO L Y,etal.Extraction algorithm for line of sight angular rate under non-Gaussiancorrelative noise environment[J].Acta Aeronautica et Astronautica Sinica,2020,41(6):324629.中的方法进行。

对比例2

进行与实施例1相同的实验，区别在于，通过以下步骤获得火箭弹视线角速率：S1、建立大长细比火箭弹弹性辨识模型；S2、利用无迹卡尔曼滤波法实时估计获得火 箭弹弹性形变量；S3、设置视线角速率提取状态模型，以弹性形变量修正导引头体视 线角，以修正后的导引头体视线角作为视线角速率提取的观测变量之一，获得视线角速率提取观测模型，基于视线角速率提取状态模型和观测模型进行滤波估计，获得视 线角速率。

在S1中，所述大长细比火箭弹弹性辨识模型基于拉格朗日函数构建，表示为：

其中，μ₁取0.05，

设置火箭弹相对地面坐标系姿态变化的运动学模型对参数θ、ψ、γ、ω_x、ω_y、 ω_z进行约束，所述火箭弹相对地面坐标系姿态变化的运动学模型表示为：

设置绕质心动力学模型对参数ω_x、ω_y、ω_z进行约束，所述绕质心动力学模型表示为：

b₁＝-σ₁η_z,1；

c₁＝-σ₁η_y,1；

a₂＝-σ₁η_y,1；

a₃＝-σ₁η_z,1

进一步地，设置瞬态坐标系下的大长细比火箭弹的动力学模型对参数M_x、M_y、 M_z、η_y、η_z、N_y、N_z、ω_x、ω_y、ω_z进行约束，所述瞬态坐标系下的大长细比火箭弹 的动力学模型表示为：

所述拉格朗日函数第1阶模态下的振动方程表示为：

/>

其中，λ₁、σ₁、M₁通过火箭弹的动能T获得。

所述火箭弹的动能T表示为：

其中，

其中截面积S＝0.07068m²，弹长l＝7.6m。

当n＝1时，通过火箭弹的动能T即可获得λ₁、σ₁、M₁。

设置全弹气动力和力矩模型对参数F_y、F_z、M_y、M_z、η_y、η_z进行约束，所述全 弹气动力和力矩模型表示为：

α(x_b,t)和β(x_b,t)可以表示为：

/>

在S2中，根据S1中获得的大长细比火箭弹弹性辨识模型，可得到耦合了弹性形变的陀螺仪角速度输出ω_xm,ω_ym,ω_zm，将陀螺仪角速度输出作为观测量 [z₁、z₂、z₃]^T＝[ω_xm,ω_ym,ω_zm]^T，获得观测模型：

利用无迹卡尔曼滤波法对观测方程进行观测，获得火箭弹弹性形变量状态的估计。

在S3中，所述视线角速率提取模型表示为：

以导引头位置处的弹性形变角对导引头体视线角进行修正，以修正后的导引头体视线角和1阶模态下的振动方程广义坐标作为视线角速率提取模型的观测变量，获得 视线角速率提取模型的观测模型，所述观测变量表示为：

[z₁,z₂,z₃,z₄]^T＝[q_αm,q_βm,η_y,η_z]^T

所述视线角速率提取模型的观测模型表示为：

采用无迹卡尔曼滤波法对观测模型进行估计，获得最终的火箭弹视线角速率。

其中，弹性辨识模型获得弹性形变量后，传递至视线角速率提取状态模型、观测模型过程中设置50ms的传输时延。

实验例1

对比实施例1与对比例1中的获取的火箭弹视线角速率，结果如图3～6所示，

图3示出了实施例1与对比例1中视线高低角估计曲线，

图4示出了实施例1与对比例1中视线方位角估计曲线，

图5示出了实施例1与对比例1中视线高低角估计误差曲线，

图6示出了实施例1与对比例1中视线方位角估计误差曲线,

从图中可以看出，实施例1能以较高的精度提取视线角速率和视线角信息，并且末端最大误差仅在0.2°以内，而对比例1中不考虑火箭弹的弹性形变而直接进行提取的方式，会导致视线角和视线角速率估计误差大幅增加，其视线角速率估计误差最大 达到了1°以上。

因此在火箭弹存在弹性形变的情况下，传统视线角速率提取方法(对比例1)不 再适用。

实验例2

对比实施例1与对比例2中获取的火箭弹视线角速率，结果如图7～9所示，其中，

图7示出了实施例1与对比例2中视线高低角估计曲线，图8示出了实施例1与对比例2中视线方位角估计曲线，图9示出了实施例1与对比例2中视线高低角估计误差曲线，

从图7～9中可以看出，实施例1与对比例2在仿真前半段弹目距离较大，弹目视线角变化较慢的情况下，都具有较高的估计精度，但是到了仿真后半段，随着弹目距离 接近，弹目视线角变化更加剧烈，对比例2由于先辨识出弹性形变量，再利用获取的弹性变形量以及将弹性变形量进行火箭弹视线角速率提取的过程具有50ms的传递时 延，导致对比例2不能像实施例1那样准确跟踪弹体的实时变化，使得对比例2中获得 的视线角速率相对于实施例有更大的误差，从而影响最终的制导精度。

在本发明的描述中，需要说明的是，术语“上”、“下”、“内”、“外”、“前”、“后” 等指示的方位或位置关系为基于本发明工作状态下的方位或位置关系，仅是为了便于 描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、 以特定的方位构造和操作，因此不能理解为对本发明的限制。此外，术语“第一”、“第二”、“第三”、“第四”仅用于描述目的，而不能理解为指示或暗示相对重要性。

在本发明的描述中，需要说明的是，除非另有明确的规定和限定，术语“安装” “相连”“连接”应作广义理解，例如，可以是固定连接，也可以是可拆卸连接，或一 体的连接普通；可以是机械连接，也可以是电连接；可以是直接连接，也可以通过中 间媒介间接连接，可以是两个元件内部的连通。对于本领域的普通技术人员而言，可以具体情况理解上述术语在本发明中的具体含义。

以上结合了优选的实施方式对本发明进行了说明，不过这些实施方式仅是范例性的，仅起到说明性的作用。在此基础上，可以对本发明进行多种替换和改进，这些均 落入本发明的保护范围内。

Claims (2)

1.一种大长细比火箭弹弹性辨识与视线角速率提取融合方法，其特征在于，包括以下步骤：

S1、以视线角、视线角速率、姿态角、角速度、火箭弹振动方程广义坐标及其广义坐标的变化率作为状态变量，建立大长细比火箭弹视线角速率状态模型；

S2、对火箭弹体视线高低角、火箭弹体视线方位角进行修正；

S3、以火箭弹体视线高低角、火箭弹体视线方位角以及火箭弹弹体三个轴向的角速度作为观测变量，建立大长细比火箭弹视线角速率观测模型；

S4、基于视线角速率状态模型和观测模型，通过滤波估计获得视线角速率；

在S1中，所述大长细比火箭弹视线角速率状态模型表示为：

其中，状态变量[x₁,x₂,x₃,x₄,x₅,x₆,x₇,x₈,x₉,x₁₀,x₁₁,x₁₂,x₁₃,x₁₄]^T设置为：

q_z表示火箭弹视线高低角，q_y表示火箭弹视线方位角；

θ、ψ、γ为地面固连坐标系下火箭弹的姿态角，θ为俯仰角，ψ为偏航角，γ为滚转角；

ω_x、ω_y、ω_z为瞬态坐标系下火箭弹的角速度；

η_y为瞬态坐标系下火箭弹振动方程拉格朗日函数第1阶模态下的垂直方向广义坐标、η_z为瞬态坐标系下火箭弹振动方程拉格朗日函数第1阶模态下的侧向广义坐标；

r表示火箭弹与目标之间的相对距离；

表示火箭弹在体视线坐标系内垂直方向加速度分量，/>

表示火箭弹在体视线坐标系内水平加速度分量；

M_x、M_y、M_z为火箭弹质心处气动弹性力矩；

J_x为火箭弹极转动惯量，J_y与J_z为赤道转动惯量；

N_y、N_z分别为瞬态坐标系下垂直方向和侧向广义力；

M₁为瞬态坐标系下火箭弹动能的拉格朗日函数第1阶模态的质量，通过火箭弹的动能T获得；

μ₁为瞬态坐标系下火箭弹势能的拉格朗日函数第1阶模态的临界阻尼系数，

ω₁为瞬态坐标系下火箭弹势能的拉格朗日函数第1阶模态的固有频率；

设置火箭弹相对地面坐标系姿态变化的运动学模型对参数θ、ψ、γ、ω_x、ω_y、ω_z进行约束，所述火箭弹相对地面坐标系姿态变化的运动学模型表示为：

其中，t表示时间；

设置绕质心动力学模型对参数ω_x、ω_y、ω_z进行约束，所述绕质心动力学模型表示为：

其中，a₁、a₂、a₃、b₁、b₂、b₃、c₁、c₂、c₃、d₁、d₂、d₃为中间变量，

b₁＝-σ₁η_z,1；

c₁＝-σ₁η_y,1；

a₂＝-σ₁η_y,1；

a₃＝-σ₁η_z,1

η_y,1表示垂直方向弹性形变的广义坐标，η_z,1表示侧向弹性形变的广义坐标，m为火箭弹质量，σ₁为常数，λ₁为常数；

设置瞬态坐标系下的大长细比火箭弹的动力学模型对参数M_x、M_y、M_z、η_y、η_z、N_y、N_z、ω_x、ω_y、ω_z进行约束，所述瞬态坐标系下的大长细比火箭弹的动力学模型表示为：

其中，V_x、V_y和V_z分别为火箭弹在地面固连坐标系下的速度分量，L为拉格朗日函数，表示为L＝T-P，T为火箭弹动能，P为火箭弹势能，D为阻尼耗散能，F_x、F_y、F_z为火箭弹所受合力在瞬态坐标系下的分量；

拉格朗日函数第1阶模态下的振动方程表示为：

所述火箭弹的动能T表示为：

其中，

M_i＝∫_mΦ_iΦ_idm，λ_i＝∫_mΦ_idm，

当n＝1时，通过火箭弹的动能T获得λ₁、σ₁、M₁；

在S2中，以导引头位置处的弹性形变角对导引头体视线角进行修正，表示为：

q_αm＝q_α-Φ′₁(x_head)η_y，q_βm＝q_β-Φ′₁(x_head)η_z，

q_αm为修正后的导引头垂直方向体视线角、q_βm为修正后的导引头水平方向体视线角，

q_α为体视线坐标系下刚性运动导引头垂直方向的体视线角、q_β为体视线坐标系下刚性运动导引头水平方向的体视线角，均由导引头输出，

Φ′₁(x_head)η_y为导引头位置处的弹性形变角的垂直分量、Φ′₁(x_head)η_z为导引头位置处的弹性形变角的侧向分量，x_head表示导引头位置，

设置全弹气动力和力矩模型对导引头位置处的弹性形变角Φ′₁(x_head)η_y、Φ′₁(x_head)η_z进行约束；

所述全弹气动力和力矩模型表示为：

其中，x_b为火箭弹轴段微元，

为火箭弹轴段微元x_b的升力系数导数、/>

为火箭弹轴段微元x_b的侧向力系数导数，Q为动压，S为火箭弹截面积，l表示火箭弹弹长，x_G为火箭弹末端到火箭弹质心的轴向距离；α(x_b,t)为在弹性形变影响下，t时刻在弹体轴向任一微元x_b的局部攻角，β(x_b,t)为在弹性形变影响下，t时刻在弹体轴向任一微元x_b的侧滑角，Φ′₁(x)表示Φ₁(x)的对x的导数；

在S3中，所述大长细比火箭弹视线角速率观测模型表示为：

其中，观测变量[z₁,z₂,z₃,z₄z₅]^T设置为：

[z₁,z₂,z₃,z₄z₅]^T＝[q_αm,q_βm,ω_xm,ω_ym,ω_zm]^T；

x_imu表示火箭弹陀螺仪安装位置，ω_xm,ω_ym,ω_zm为火箭弹耦合弹性形变后的陀螺仪角速度的输出，a、b、c为中间量，

为地面坐标系到瞬态坐标系的转换矩阵。

2.根据权利要求1所述的大长细比火箭弹弹性辨识与视线角速率提取融合方法，其特征在于，

S4中，基于视线角速率状态模型和观测模型，通过无迹卡尔曼滤波法进行滤波估计，获得视线角速率。

Images