基于SVD-CKF的战术导弹气动参数在线辨识方法研究

Research on Tactical Missile Aerodynamic Parameter Online Identification Method Based on SVD-CKF

摘要

在导弹气动参数辨识领域，传统扩展卡尔曼滤波（extended Kalman filter，EKF）算法往往计算复杂、计算精度低，且求解系统雅各比矩阵难.针对这个问题，本文提出了一种基于SVD-CKF的战术导弹气动参数在线辨识方法.利用容积卡尔曼滤波（cubature Kalman filter，CKF）的容积点线性化特性，避免了对雅各比矩阵的直接求解，从而降低了计算复杂度.同时，通过引入奇异值分解（singular value decomposition，SVD）技术，有效解决了传统CKF算法中可能导致协方差矩阵负定的情况，进一步提升了滤波稳定性.仿真结果表明，在六自由度战术导弹气动参数在线辨识问题中，SVD-CKF算法展现更高的辨识精度、更快的收敛速度以及更强的鲁棒性.

Abstract

In the field of missile aerodynamic parameter identification， traditional extended kalman filter （EKF） algorithms often encounter issues such as high computational complexity， low accuracy， and difficulties in solving the system’s Jacobian matrix. To address these challenges， an online identification method for missile aerodynamic parameters based on singular value decomposition-cubature kalman filter （SVD-CKF） is proposed. Leveraging the cubature point linearization characteristic of CKF， this method avoids the direct solution of the Jacobian matrix， thereby reducing computational complexity. Additionally， by introducing Singular Value Decomposition （SVD） technology， it effectively resolves the issue of potential negative definiteness in the covariance matrix that may arise in traditional CKF algorithms， further enhancing filter stability. Simulation results demonstrate that in the context of online identification of aerodynamic parameters for six-degree-of-freedom tactical missiles， the SVD-CKF algorithm exhibits higher identification accuracy， faster convergence speed， and stronger robustness.

Graphical abstract

关键词

在线参数辨识 / 容积卡尔曼滤波 / 奇异值分解 / 导弹

Key words

on-line parameter identification / cubature Kalman filtering / singular value decomposition / missile

导弹研制过程中，其动态模型的构建主要依赖于三种方法：理论计算^［1］、风洞试验^［2］以及飞行试验^［3］.在导弹设计初期，利用理论计算来预测导弹的动态行为，但由于空气动力学理论尚不完全成熟，计算结果往往与实际情况存在显著偏差.在导弹设计的中后期，通过风洞试验在实验室环境中模拟导弹所受的空气动力学，以验证和修正理论模型.由于通常采用缩比模型模拟实际导弹，故难以完全复现真实飞行条件.相较而言，在导弹设计的最后阶段，进行飞行试验能够更贴近实际飞行场景，通过在实际飞行中观察和记录导弹的实际行为来验证理论和风洞试验的结果，但其成本高昂且不易重复进行.鉴于理论计算存在的偏差以及风洞试验和飞行试验的局限性，采用参数辨识获得导弹的气动特性对理论计算进行修正，验证风洞试验的结果，并优化飞行试验的设计，此举可为后续的高精度自适应控制奠定坚实的基础^［4］.

在参数辨识中有离线辨识和在线辨识，如最小二乘法、递推最小二乘法、带遗忘因子的最小二乘法等.在线辨识常常使用卡尔曼滤波^［5］、扩展卡尔曼滤波、无迹卡尔曼滤波等.相对于离线辨识，在线辨识可以在线实时估计模型参数，用于控制系统模型或参数在线修订并进一步提高精度，故比离线辨识应用得更加普遍.

在导弹在线参数辨识领域，EKF是目前应用最为广泛的参数辨识算法，如文献［6］对导弹助推段气动系数建模与在线辨识方法进行研究，利用EKF对战术导弹进行参数辨识；但是对于复杂非线性系统，求解系统的雅各比矩阵和海森矩阵很困难，导致其算法存在计算复杂、估计精度不高的问题.为了解决估计精度不高的问题，文献［7］提出一种基于迭代滤波理论的扩展卡尔曼滤波改进算法，通过在EKF估计点上重新进行泰勒级数展开并迭代，有效减少了线性化误差，从而提高了辨识精度，但是没有从根源上解决计算复杂的问题.为解决算法计算复杂的问题，文献［8］基于无迹卡尔曼滤波（unscented Kalman filter， UKF）的巡飞导弹气动参数在线辨识方法，采用无迹卡尔曼滤波算法实现在线参数辨识，验证了其精度比EKF更高，同时有效避免了求解雅各比矩阵的问题，极大减少了计算量.

综上所述，为解决卡尔曼滤波在估计气动参数时精度低、计算复杂的问题^［9］，本文提出了一种基于SVD-CKF的战术导弹气动参数在线辨识方法.相较于EKF，SVD-CKF在处理高维度且高度非线性的系统滤波问题时，展现出更低的计算负担，且无须进行雅各比矩阵的求解^［10］.相较于传统CKF滤波算法，本文利用类似文献［11］的思想，采用数值稳定性较强的SVD分解代替传统CKF中的Cholesky分解，对协方差矩阵进行分解，可进一步提高滤波的稳定性与算法的鲁棒性^［12］.最后仿真结果表明，SVD-CKF算法在面临高噪声和复杂动态环境时，对导弹气动参数的估计精度和滤波性能较高，为导弹制导与控制系统的设计提供了更为可靠的数据支持.

1 战术导弹气动参数辨识模型

1.1 状态方程

以某型战术导弹为对象，在建模过程中，为了简化分析并突出主要影响因素，忽略地球自转和曲率对导弹运动产生的微小影响，同时假设导弹推进系统产生的推力在此阶段可忽略不计.此外，考虑到导弹的对称性，将其视为一个转动惯量积为零的系统.在弹体坐标系上建立了该战术导弹的六自由度状态方程如下所示：

{\dot{V}}_x={\omega }_z{V}_y-{\omega }_y{V}_z-g\mathrm{s}\mathrm{i}\mathrm{n}\phi -X/m

{\dot{V}}_y={\omega }_x{V}_z-{\omega }_z{V}_x-g\mathrm{c}\mathrm{o}\mathrm{s}\phi \mathrm{c}\mathrm{o}\mathrm{s}\gamma +Y/m

{\dot{V}}_z={\omega }_y{V}_x-{\omega }_x{V}_y+g\mathrm{c}\mathrm{o}\mathrm{s}\phi \mathrm{s}\mathrm{i}\mathrm{n}\gamma +Z/m

\dot{\phi }={\omega }_y\mathrm{s}\mathrm{i}\mathrm{n}\gamma +{\omega }_z\mathrm{c}\mathrm{o}\mathrm{s}\gamma

\dot{\psi }=({\omega }_y\mathrm{c}\mathrm{o}\mathrm{s}\gamma -{\omega }_z\mathrm{s}\mathrm{i}\mathrm{n}\gamma )\frac{\mathrm{1}}{\mathrm{c}\mathrm{o}\mathrm{s}\phi }

\dot{\gamma }={\omega }_x-({\omega }_y\mathrm{c}\mathrm{o}\mathrm{s}\gamma -{\omega }_z\mathrm{s}\mathrm{i}\mathrm{n}\gamma )\mathrm{t}\mathrm{a}\mathrm{n}\phi

\dot{h}=V_x\mathrm{s}\mathrm{i}\mathrm{n}\phi +V_y\mathrm{c}\mathrm{o}\mathrm{s}\phi \mathrm{c}\mathrm{o}\mathrm{s}\gamma -V_z\mathrm{c}\mathrm{o}\mathrm{s}\phi \mathrm{s}\mathrm{i}\mathrm{n}\gamma

{\dot{\omega }}_x=\frac{I_y-I_z}{I_x}{\omega }_y{\omega }_z+\frac{M_x}{I_x}

{\dot{\omega }}_y=\frac{I_z-I_x}{I_y}{\omega }_x{\omega }_z+\frac{M_y}{I_y}

{\dot{\omega }}_z=\frac{I_x-I_y}{I_z}{\omega }_x{\omega }_y+\frac{M_z}{I_z}（1）

式中：X\mathrm{、}Y\mathrm{、}Z分别为轴向力、法向力、横向力；V_x\mathrm{、}{V}_y\mathrm{、}{V}_z分别是弹体坐标系轴向、法向和侧向三个分量的线速度；\phi \mathrm{、}\gamma \mathrm{、}\psi分别为俯仰角、滚转角、偏航角；{\omega }_x\mathrm{、}{\omega }_y\mathrm{、}{\omega }_z分别是三个方向的角速度；I_x\mathrm{、}{I}_{\mathrm{y}}\mathrm{、}{I}_z分别是滚转、偏航和俯仰通道的转动惯量；M_x\mathrm{、}{M}_y\mathrm{、}{M}_z分别是所有的外力对导弹质心的力矩在弹体坐标系各轴上的分量；g为重力加速度，本文取9.8 \mathrm{m}/{\mathrm{s}}^{\mathrm{2}}；h为高度信息； m是导弹质量.

1.2 量测方程

量测量包括角速度信息、过载信息、姿态角信息、高度以及攻角和侧滑角.量测方程为

{\omega }_x={\omega }_x+v_{\mathrm{1}},{\omega }_y={\omega }_y+v_{\mathrm{2}},{\omega }_z={\omega }_z+v_{\mathrm{3}}

N_x=\frac{-X}{mg}+v_{\mathrm{4}},N_y=\frac{-Y}{mg}+v_{\mathrm{5}},N_z=\frac{Z}{mg}+v_{\mathrm{6}}

\phi =\phi +v_{\mathrm{7}},\psi =\psi +v_{\mathrm{8}},\gamma =\gamma +v_{\mathrm{9}},h=h+v_{\mathrm{10}}

\alpha =\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{t}\mathrm{a}\mathrm{n}(-V_y/V_x)+v_{\mathrm{11}},\beta =\mathrm{a}\mathrm{r}\mathrm{c}\mathrm{s}\mathrm{i}\mathrm{n}(V_z/V)+v_{\mathrm{12}}（2）

式中：\alpha \mathrm{、}\beta分别为攻角和侧滑角；N_x\mathrm{、}{N}_y\mathrm{、}{N}_z分别为各个方向上的过载；V为速度；v_{\mathrm{1}}~v_{\mathrm{12}}为服从高斯分布的观测噪声.

1.3 气动系数模型

在六自由度战术导弹系统设计中，气动参数的准确建模对导弹的飞行性能至关重要.通常，采用代数多项式构建气动参数模型，以全面捕捉导弹在不同飞行条件下的气动特性.根据参考文献［6］中的模型，本文采用了以下气动参数模型：

\begin{array}{l}{C}_x=C_{x\mathrm{0}}+C_x^{\alpha }\alpha +C_x^{\beta }\beta +C_x^{{\delta }_x}{\delta }_x\\ C_y=C_{y\mathrm{0}}+C_y^{\alpha }\alpha +C_y^{{\delta }_y}{\delta }_y\\ C_z=C_{z\mathrm{0}}+C_z^{\beta }\beta +C_z^{{\delta }_z}{\delta }_z\end{array}（3）

式中：{\delta }_x\mathrm{、}{\delta }_y\mathrm{、}{\delta }_z分别为俯仰舵、偏航舵、滚转舵；C_{x\mathrm{0}}\mathrm{、}{C}_{y\mathrm{0}}\mathrm{、}{C}_{z\mathrm{0}}为零攻角、零侧滑角下的基准气动系数；C_x^{\alpha }\mathrm{、}{C}_y^{\alpha }分别为轴向力、法向力对攻角相关的修正项；C_x^{\beta }\mathrm{、}{C}_z^{\beta }分别为法向力、横向力对侧滑角相关的修正项；C_x^{{\delta }_x}\mathrm{、}{C}_y^{{\delta }_y}\mathrm{、}{C}_z^{{\delta }_z}为各方向对控制舵偏相关的修正项.

由式（3）可知，需要辨识的气动参数为C_{x\mathrm{0}},C_x^{\alpha },C_x^{\beta },C_x^{{\delta }_x},C_{y\mathrm{0}},C_y^{\alpha },C_y^{{\delta }_y},C_{z\mathrm{0}},C_z^{\beta },C_z^{{\delta }_z}共10个，这些参数构成待辨识参数向量\mathit{\theta }为：

\mathit{\theta }={\left[C_{x\mathrm{0}},C_x^{\alpha },C_x^{\beta },C_x^{{\delta }_x},C_{y\mathrm{0}},C_y^{\alpha },C_y^{{\delta }_y},C_{z\mathrm{0}},C_z^{\beta },C_z^{{\delta }_z}\right]}^{\mathrm{T}}（4）

2 基于SVD-CKF的气动参数在线辨识方法

2.1 在线辨识过程

参数辨识过程如图1所示，基于构建的状态方程，将待辨识参数向量\mathit{\theta }纳入状态向量中进行初始化.通过状态方程和量测方程分别完成时间更新和量测更新步骤，以获取状态预测值和量测预测值.最终，得到气动参数的精确辨识结果.

2.2 SVD-CKF滤波辨识器设计

CKF的核心特色在于其采用确定性数值积分策略，摒弃了传统卡尔曼滤波在处理非线性函数时常用的线性化方法.相较于EKF通过雅可比矩阵实现非线性函数的线性近似，CKF则利用立方点技术，更为精确地评估高斯分布在传播与更新过程中的动态变化.这些立方点在状态空间中执行采样操作，并在预测与更新阶段被有效利用，以捕捉非线性系统内部的更深层次信息^［13］.

考虑下列非线性估计的离散高斯模型：

\dot{\mathit{x}}(\mathit{k})=f(\mathit{x}(k-\mathrm{1}),u(k-\mathrm{1}))+w(k-\mathrm{1})

\mathit{z}(\mathit{k})=h(\mathit{x}(k),u(k))+v(k)（5）

式中：\mathit{x}(k)为状态矢量；\mathit{z}(k)为观测矢量；u(k)为控制输入；w(k)为系统噪声；v(k)为观测噪声.这些噪声均为零均值高斯白噪声.

步骤1 滤波器初始化

设置状态矢量的初始值{\widehat{\mathit{x}}}_{\mathrm{0}}及相应的误差协方差矩阵{\mathit{P}}_{\mathrm{0}}.

\begin{array}{l}{\widehat{\mathit{x}}}_{\mathrm{0}}=E[{\mathit{x}}_{\mathrm{0}}]\\ {\mathit{P}}_{\mathrm{0}}=E[({\mathit{x}}_{\mathrm{0}}-{\widehat{\mathit{x}}}_{\mathrm{0}})({\mathit{x}}_{\mathrm{0}}-{\widehat{\mathit{x}}}_{\mathrm{0}}{)}^{\mathrm{T}}]\end{array}（6）

步骤2 计算容积点与相应的权值，利用3阶径向规则计算标准高斯加权积分，即：

E[\mathit{x}|\mathit{z}]=\underset{R}{\int }f(\mathit{x})N(\mathit{x};\mathrm{0},\mathit{I})\mathrm{d}x\approx {\displaystyle \sum _{i=\mathrm{1}}^m}f({\mathit{\xi }}_i)

{\mathit{\xi }}_i=\sqrt[]{\frac{m}{\mathrm{2}}}{[\mathrm{1}]}_i,i=\mathrm{1,2},\cdots ,m=\mathrm{2}n

{\omega }_i=\frac{\mathrm{1}}{m},i=\mathrm{1,2},\cdots ,m=\mathrm{2}n（7）

式中：m代表容积采样点个数，使用3阶容积规则时容积采样点是状态向量维数的2倍；［1］代表n维空间的点集，即：

\left\{\left(\begin{array}{l}\mathrm{1}\\ \mathrm{0}\\ ⋮\\ \mathrm{0}\end{array}\right)\cdots \left(\begin{array}{l}\mathrm{0}\\ \mathrm{0}\\ ⋮\\ \mathrm{1}\end{array}\right)\left(\begin{array}{l}-\mathrm{1}\\ \mathrm{0}\\ ⋮\\ \mathrm{0}\end{array}\right)\cdots \left(\begin{array}{l}\mathrm{0}\\ \mathrm{0}\\ ⋮\\ -\mathrm{1}\end{array}\right)\right\}（8）

步骤3 时间更新

在此步骤中，若采用传统的Cholesky分解，CKF算法进行多次滤波之后，容易导致协方差矩阵出现负定的情况，使得滤波发散或截止^［14］.本文提出使用SVD分解代替传统Cholesky分解，如式（9）所示：

{\mathit{P}}_{k-\mathrm{1}|k-\mathrm{1}}=\mathit{U}\left(\begin{array}{cc}\mathit{S}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}\end{array}\right){\mathit{V}}^{\mathrm{T}}（9）

式中，\mathit{S}为对角矩阵，经过协方差矩阵的SVD分解：

{\mathit{P}}_{k-\mathrm{1}|k-\mathrm{1}}={\mathit{U}}_{k-\mathrm{1}}{\mathit{S}}_{k-\mathrm{1}}{\mathit{V}}_{k-\mathrm{1}}^{\mathrm{T}}

{\mathit{X}}_{i,k|k-\mathrm{1}}={\mathit{U}}_{i,k-\mathrm{1}}\sqrt[]{{\mathit{S}}_i}{\xi }_i+{\widehat{\mathit{x}}}_{k|k-\mathrm{1}}（10）

计算通过非线性状态方程传播Cubature点：

{\mathit{X}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{*}}=f({\mathit{X}}_{i,k|k-\mathrm{1}},{\mathit{u}}_{k-\mathrm{1}})（11）

计算状态预测值和状态误差协方差矩阵的值：

{\widehat{\mathit{x}}}_{i,k|k-\mathrm{1}}=\frac{\mathrm{1}}{\mathrm{2}n}{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{2}n}}{\mathit{X}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{*}}

{\mathit{P}}_{k|k-\mathrm{1}}=\frac{\mathrm{1}}{\mathrm{2}n}{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{2}n}}{\mathit{X}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{*}}{\mathit{X}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{*}\mathrm{T}}-{\widehat{\mathit{x}}}_{k|k-\mathrm{1}}{{\widehat{\mathit{x}}}^{\mathrm{T}}}_{k|k-\mathrm{1}}+{\mathit{Q}}_{k-\mathrm{1}}（12）

步骤4 量测更新

对{\mathit{P}}_{k|k-\mathrm{1}}进行SVD分解，并计算容积点：

{\mathit{P}}_{k|k-\mathrm{1}}={\mathit{U}}_{k|k-\mathrm{1}}{\mathit{S}}_{k|k-\mathrm{1}}{\mathit{V}}_{{}_{k|k-\mathrm{1}}}^{\mathrm{T}}

{\mathbf{X}}_{i,k|k-\mathrm{1}}={\mathit{U}}_{i,k|k-\mathrm{1}}\sqrt[]{{\mathit{S}}_i}{\mathit{\xi }}_i+{\widehat{\mathit{x}}}_{k|k-\mathrm{1}}（13）

通过非线性量测方程计算传播的容积点：

{\mathit{Z}}_{i,k|k-\mathrm{1}}=h({\mathit{X}}_{i,k|k-\mathrm{1}})（14）

对非线性量测预测值、新协方差和量测协方差矩阵进行计算：

\left\{\begin{array}{l}{\widehat{\mathit{Z}}}_{i,k|k-\mathrm{1}}=\frac{\mathrm{1}}{\mathrm{2}n}{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{2}n}}{\mathit{Z}}_{i,k|k-\mathrm{1}}\\ {\mathit{P}}_{zz,k|k-\mathrm{1}}=\frac{\mathrm{1}}{\mathrm{2}n}{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{2}n}}{\mathit{Z}}_{i,k|k-\mathrm{1}}{\mathit{Z}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{T}}-{\widehat{\mathit{z}}}_{k|k-\mathrm{1}}{{\widehat{\mathit{z}}}^{\mathrm{T}}}_{k|k-\mathrm{1}}+{\mathit{R}}_k\\ {\mathit{P}}_{xz,k|k-\mathrm{1}}=\frac{\mathrm{1}}{\mathrm{2}n}{\displaystyle \sum _{i=\mathrm{1}}^{\mathrm{2}n}}{\mathit{X}}_{i,k|k-\mathrm{1}}{\mathit{Z}}_{{}_{i,k|k-\mathrm{1}}}^{\mathrm{T}}-{\widehat{\mathit{x}}}_{k|k-\mathrm{1}}{{\widehat{\mathit{z}}}^{\mathrm{T}}}_{k|k-\mathrm{1}}\end{array}\right.（15）

计算滤波增益矩阵、状态与误差协方差矩阵：

{\mathit{W}}_k={\mathit{P}}_{xz,k|k-\mathrm{1}}{\mathit{P}}_{xz,k|k-\mathrm{1}}^{-\mathrm{1}}

{\widehat{\mathit{x}}}_{k|k}={\widehat{\mathit{x}}}_{k|k-\mathrm{1}}+{\mathit{W}}_k({\mathit{z}}_k-{\widehat{\mathit{z}}}_{k|k-\mathrm{1}})

{\mathit{P}}_{k|k}={\mathit{P}}_{k|k-\mathrm{1}}-{\mathit{W}}_k{\mathit{P}}_{zz,k|k-\mathrm{1}}{\mathit{W}}_k^{\mathrm{T}}（16）

3 仿真分析

卡尔曼滤波是基于贝叶斯估计的线性动态系统状态估计算法，其核心在于实现状态的最优估计.在面临未知参数辨识的复杂问题时，单纯的卡尔曼滤波算法显然无法直接应用.为此，将待辨识的参数增广到状态向量中，状态向量扩展如式（17）所示：

\mathit{x}=[V_x,V_y,V_z,{\omega }_x,{\omega }_y,{\omega }_z,\phi ,\psi ,\gamma ,h,C_{x\mathrm{0}},C_x^{\alpha },C_x^{\beta },

C_x^{{\delta }_x},C_{y\mathrm{0}},C_y^{\alpha },C_y^{{\delta }_y},C_{z\mathrm{0}},C_z^{\beta },C_z^{{\delta }_z}]（17）

3.1 仿真条件

选择某助推滑翔战术导弹进行研究.导弹的质量m=\mathrm{980}\mathrm{k}\mathrm{g}，特征参考面积S=\mathrm{1.677}\mathrm{1}{\mathrm{m}}^{\mathrm{2}}，特征长度L=\mathrm{3.85}\mathrm{m}，导弹的标准气动数据是攻角与马赫数的函数，每一步的计算周期T=\mathrm{0.005}\mathrm{s}.取初始速度为{(\mathrm{1}\mathrm{000,10,10})}^{\mathrm{T}}\mathrm{m}/\mathrm{s}，初始角速度为{(\mathrm{3}°,\mathrm{3}°,\mathrm{3}°)}^{\mathrm{T}}/\mathrm{s}，初始姿态角为{(\mathrm{15}°,\mathrm{3}°,\mathrm{3}°)}^{\mathrm{T}}/\mathrm{s}.

设定式（2）量测方程中过载的测量误差为{\mathrm{10}}^{-\mathrm{5}}\mathrm{m}/{\mathrm{s}}^{\mathrm{2}}，姿态角的测量误差为{\mathrm{10}}^{-\mathrm{5}}(°)/\mathrm{s}，角速度的测量误差为{\mathrm{10}}^{-\mathrm{5}}\mathrm{r}\mathrm{a}\mathrm{d}/\mathrm{s}，高度的测量误差为{\mathrm{10}}^{-\mathrm{2}}\mathrm{m}/\mathrm{s}.量测噪声误差符合高精度仿真环境的实验要求.

根据实际传感器的测量误差协方差矩阵，设定量测噪声的协方差矩阵是一个12×12的对角矩阵，如式（18）所示：

\mathit{R}=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}[\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0.000}\mathrm{1};\mathrm{0.000}\mathrm{1};\mathrm{0.000}\mathrm{1};

\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{1}]\times {\mathrm{10}}^{-\mathrm{3}}（18）

系统过程噪声协方差矩阵为20×20的对角矩阵

Q=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}\left[\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1};\mathrm{1}\right]\times {\mathrm{10}}^{-\mathrm{5}}（19）

初始的协方差矩阵为20×20的对角矩阵，具体如式（20）所示：

P=\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}[\mathrm{4};\mathrm{0.000}\mathrm{4};\mathrm{0.000}\mathrm{4};\mathrm{0};\mathrm{0};\mathrm{0};

\mathrm{0.006}\mathrm{4};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};\mathrm{0};

\mathrm{0.000}\mathrm{4};\mathrm{0};\mathrm{0};\mathrm{0.000}\mathrm{3};\mathrm{0}]\times {\mathrm{10}}^{\mathrm{4}}（20）

3.2 输入信号设计

为了更好地了解导弹的空气动力学特性，需要设定特定的控制输入信号来激发系统的不同运动模态，从而获得足够的信息进行参数辨识^［15］.因此输入信号应尽量包含多种频率成分，满足要求的常用输入信号有阶跃输入、方波输入、偶极方波输入等.各个方向舵偏角输入激励如图2所示.

3.3 辨识结果

本文精选了如图3~图8所示的典型参数估计结果.SVD-CKF迭代中避免了先验协方差矩阵负定导致的滤波发散，确保了滤波器运行.相较于CKF，SVD-CKF收敛速度更快，滤波更新更稳定.由图5发现有些参数的辨识结果在真值附近波动，存在一定的辨识误差，这可能与激励信号的设计有关.对比分析表1和表2的辨识结果，可以得出结论：与CKF算法相比，SVD-CKF在辨识精度上更优.具体而言，CKF算法在气动参数辨识中的相对误差保持在5%以内，而SVD-CKF则将所有辨识结果的相对误差控制在3%以内，显现出更优的辨识性能.

4 结 论

为解决EKF计算复杂且精度低的问题，提出了基于SVD-CKF的战术导弹气动参数在线辨识方法.通过采样点直接近似非线性函数，减少了计算量；提出使用SVD分解替代Cholesky分解，增强滤波稳定性.仿真结果表明，相较于CKF算法，SVD-CKF能够将气动参数辨识的相对误差控制在3%以内，提升了2个百分点的辨识精度，气动参数辨识精度得到优化；SVD-CKF算法还表现出了更快的收敛速度和更高的滤波稳定性.

参考文献

[1] TISCHLER M B， REMPLE R K. Aircraft and rotorcraft system identification： engineering methods and flight test examples［M］. Beijing： Aviation Industry Press， 2012.

[2] 邹新生 ．飞行器非线性参数辨识与鲁棒控制研究［D］．北京：清华大学，2006．

[3] ZOU X S. Research on nonlinear parameter identification and robust control of aircraft［D］．Beijing：Tsinghua University，2006．（in Chinese）

[4] 张淼，刘铁军，马涂亮，等 ．基于CFD方法的大型客机高速气动设计［J］．航空学报，2016，37（1）： 244-254．

[5] ZHANG M， LIU T J， MA T L，et al ．High speed aerodynamic design of large civil transporter based on CFD method［J］．Acta Aeronautica et Astronautica Sinica，2016，37（1）：244-254．（in Chinese）

[6] 刘斌，刘沛清，王亮 ．某小型边条翼无人飞机的气动参数估算及风洞试验分析［J］．飞机设计，2010，30（3）：1-5．

[7] LIU B， LIU P Q， WANG L ．A small trim-wing unmanned aircraft aerodynamic parameter estimation and analysis of wind tunnel test［J］．Aircraft Design，2010，30（3）：1-5．（in Chinese）

[8] KALMAN R E ．A new approach to linear filtering and prediction problems［J］．Journal of Basic Engineering，1960，82（1）：35-45．

[9] 孙乐园， 江振宇， 黄振新. 导弹助推段气动系数建模与在线辨识方法研究［C］//2016 IEEE中国制导、导航与控制学术会议.2016：2574-2579.

[10] SUN L Y， JIANG Z Y， HUANG Z X. Research on modeling and online identification method of aerodynamic coefficients for missile booster segments［C］//Proceedings of the 2016 IEEE Chinese Guidance， Navigation and Control Conference. 2016：2574-2579. （in Chinese）.

[11] 罗鹏，杨华，陈伟芳 ．一种改进的扩展卡尔曼滤波气动参数辨识方法［J］．工业控制计算机，2018，31（12）：13-16．

[12] LUO P， YANG H， CHEN W F ．An improved method of extended Kalman filter for aerodynamic parameter identification［J］．Industrial Control Computer，2018，31（12）：13-16．（in Chinese）

[13] 涂海峰，贾生伟，阳丰俊，等 ．基于无迹卡尔曼滤波的巡飞弹气动参数在线辨识方法［J］．航天控制，2018，36（5）：14-18．

[14] TU H F， JIA S W， YANG F J，et al ．Loitering munition aerodynamic parameter online identification method based on unscented Kalman filter［J］．Aerospace Control，2018，36（5）：14-18．（in Chinese）

[15] MALLICK M， TIAN X Q， LIU J ．Evaluation of measurement converted KF，EKF，UKF，CKF，and PF in GMTI filtering［C］//2021 International Conference on Control，Automation and Information Sciences （ICCAIS）． Xi’an，China． IEEE，2021：21-27．

[16] 刘世林， 李德俊， 姚伟， 等 ．基于核极限学习机与容积卡尔曼滤波融合的锂电池荷电状态估计［J］．湖南大学学报（自然科学版）， 2023， 50（10）： 51-59．

[17] LIU S L， LI D J， YAO W，et al ．Estimation on state of charge of lithium battery based on fusion of kernel extreme learning machine and cubature Kalman filter［J］．Journal of Hunan University （Natural Sciences）， 2023， 50（10）： 51-59．（in Chinese）

[18] 王姚宇， 陈仁文， 张祥 ．改进的基于奇异值分解的抗差容积卡尔曼滤波算法在全球定位导航中的应用［J］．科学技术与工程，2021，21（6）：2356-2362．

[19] WANG Y Y， CHEN R W， ZHANG X ．Application of improved robust singular value decomposition-cubature Kalman filter algorithm in global positioning system navigation［J］．Science Technology and Engineering，2021，21（6）：2356-2362．（in Chinese）

[20] WANG D， YANG J， XIONG K，et al ．CNS/INS integrated navigation method based on improved adaptive CKF algorithm［C］//IECON 2023- 49th Annual Conference of the IEEE Industrial Electronics Society． Singapore，Singapore． IEEE，2023：1-7．

[21] VERMA H O， PEYADA N K ．Estimation of aerodynamic parameters near stall using maximum likelihood and extreme learning machine-based methods［J］．The Aeronautical Journal，2021，125（1285）：489-509．

[22] ZHAO J F， ZHAO X， LI Q H，et al ．Research on ultra-tight integrated navigation method based on CKF filtering algorithm［J］．Journal of Physics：Conference Series，2022，2221（1）：012048．

[23] WEI C Z， LV J X， LI Y L，et al ．SVR-ND method for online aerodynamic parameter estimation［J］．IEEE Access，2020，8：207204-207215．

基金资助

国家自然科学基金资助项目(62473138)

National Natural ScienceFoundation of China(62473138)

广东省自然科学基金青年提升计划项目(2024A1515030184)

广州市增城区重点研发项目(2024ZCKJ01)

湖南省自然科学基金资助项目(2022JJ30162)

Natural Science Foundation of Hunan Province(2022JJ30162)