CN111766397A - 一种基于惯性/卫星/大气组合的气象风测量方法 - Google Patents

一种基于惯性/卫星/大气组合的气象风测量方法 Download PDF

Info

Publication number
CN111766397A
CN111766397A CN202010563456.4A CN202010563456A CN111766397A CN 111766397 A CN111766397 A CN 111766397A CN 202010563456 A CN202010563456 A CN 202010563456A CN 111766397 A CN111766397 A CN 111766397A
Authority
CN
China
Prior art keywords
wind
speed
axis
angle
aircraft
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010563456.4A
Other languages
English (en)
Other versions
CN111766397B (zh
Inventor
王学运
从浩
张京娟
赵明
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hefei Innovation Research Institute of Beihang University
Original Assignee
Hefei Innovation Research Institute of Beihang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hefei Innovation Research Institute of Beihang University filed Critical Hefei Innovation Research Institute of Beihang University
Priority to CN202010563456.4A priority Critical patent/CN111766397B/zh
Publication of CN111766397A publication Critical patent/CN111766397A/zh
Application granted granted Critical
Publication of CN111766397B publication Critical patent/CN111766397B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P5/00—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01D—MEASURING NOT SPECIALLY ADAPTED FOR A SPECIFIC VARIABLE; ARRANGEMENTS FOR MEASURING TWO OR MORE VARIABLES NOT COVERED IN A SINGLE OTHER SUBCLASS; TARIFF METERING APPARATUS; MEASURING OR TESTING NOT OTHERWISE PROVIDED FOR
    • G01D21/00—Measuring or testing not otherwise provided for
    • G01D21/02—Measuring two or more variables by means not covered by a single other subclass
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P13/00—Indicating or recording presence, absence, or direction, of movement
    • G01P13/02—Indicating direction only, e.g. by weather vane
    • G01P13/025—Indicating direction only, e.g. by weather vane indicating air data, i.e. flight variables of an aircraft, e.g. angle of attack, side slip, shear, yaw

Abstract

本发明公开一种基于惯性/卫星/大气组合的气象风测量方法，尤其涉及风速风向的在线测量，在不改变飞行器结构外形、不增加测量装置和硬件设备的前提下，利用现有的机载惯性导航系统、卫星定位系统、大气数据系统的输出参数，得到惯性系统的加速度、角速度和姿态角，卫星测得地速、大气数据系统测得真空速，通过构建无迹卡尔曼滤波器，将风速风向测量模型和状态估计过程融合起来，实现飞行器机体周围风速风向的实时测量，并完成下一时刻风速风向的估计。本发明所提出的方法方案可以解决飞行器周围局域风速风向的实时测量和大区域连续飞行实时风速风向测量问题，为四维航迹飞行、飞机准点达到提供方法支持。

Description

一种基于惯性/卫星/大气组合的气象风测量方法
方法领域
本发明涉及风速风向测量方法领域，尤其涉及一种基于惯性/卫星/大气组合的气象风测量方法，可以解决飞行器周围局域风速风向的实时测量和大区域连续飞行实时风速风向测量问题，为四维航迹飞行、飞机准点达到提供方法支持。
背景方法
飞行器的飞行状态不可避免的受机体周围风速风向的影响，风速风向的精确测量对提高飞行的飞行可靠性、安全性和飞行品质意义重大。目前主要的测风方法是单一的机载大气数据系统提供是风速风向信息，或大气数据系统与卫星定位系统进行组合测量。其中大气数据系统由于测量传感设备存在测量迟滞，所以单一的大气数据系统测量不能实现机动情况下实时测量；卫星定位系统测量时，卫星信号易受干扰，还可能会发生信号中断等情况，都不能满足飞行器周围局域风速风向的实时测量和大区域连续飞行实时风速风向的测量问题。
发明内容
本发明目的就是为了弥补已有方法的缺陷，提供一种基于惯性/卫星/大气组合的气象风测量方法。
本发明是通过以下方法方案实现的：
一种基于惯性/卫星/大气组合的气象风测量方法，利用飞行器机载航电设备输出的参数，包括：
惯性系统输出的加速度、角速度和姿态角；
卫星定位系统输出的地速；
大气数据系统输出的真空速；
由惯性导航系统的输出信息，对状态转换矩阵进行计算；
由卫星定位系统输出的地速通过机体坐标系系和地理坐标系的状态转移矩阵，求解出机体坐标系下的飞行器三轴速度。
由大气数据系统测量数据，对飞行器的真空速进行解算。
由当前时刻飞行器的地速，状态转移矩阵和真空速解算结果，对飞行器机身周围的风速风向进行计算，并通过无迹卡尔曼滤波器对下一时刻风速风向进行实时精确估计。
包含以下步骤：
第一步，以周期△T读取惯性导航系统输出的三个角速度、三个姿态角、三个线加速度信息，三个角速度[p q r]T分别是机体系下X轴、Y轴和Z轴方向上的横滚角速度p，俯仰角速度q，偏航角速度r；三个姿态角[φ θ ψ]T分别为横滚角φ、俯仰角θ和偏航角ψ；三个线加速度[ax ay az]T分别为机体坐标系下X轴方向的线加速度ax、Y轴方向的线加速度ay和Z轴方向的线加速度az，其中机体坐标系的X轴、Y轴、Z轴的指向分别为向前、向右和向下。
第二步，以周期△T读取机载卫星定位系统输出的地理坐标系下的三轴地速[Vx VyVz]T，其中Vx为飞行器在地理坐标系X轴的速度分量、Vy为飞行器在地理坐标系Y轴的速度分量、Vz为飞行器在地理坐标系的速度分量，其中地理坐标系的X轴、Y轴、Z轴的指向分别为北向、东向和地。
第三步，状态转移矩阵为地理坐标系到机体坐标系的状态转换矩阵，通过姿态角信息获得，姿态角信息包括横滚角φ、俯仰角θ和偏航角ψ，计算状态转移矩阵
Figure BDA0002546877930000021
Figure BDA0002546877930000022
第四步，通过状态转移矩阵可以地理坐标系下的三轴地速[Vx Vy Vz]T转换到机体系下：
Figure BDA0002546877930000023
[wx wy wz]T为三轴风速。
第五步，机载真空速通常由机载大气数据系统中压力传感器及温度传感器提供的压力以及温度，结合伯努利方程解算得到的，计算公式：
Figure BDA0002546877930000024
其中，△P为动压，ρ为空气密度，K为校正因子。
第六步，利用卫星定位系统测得的地速
Figure BDA0002546877930000025
惯性导航系统测得的姿态角，大气数据系统提供的攻角α和侧滑角β解算得到的转态转移矩阵、大气数据系统测得的真空速
Figure BDA0002546877930000031
得到：
Figure BDA0002546877930000032
其中
Figure BDA0002546877930000033
wx为X轴方向上的风速，wy为Y轴方向上的风速，wz是Z轴方向上的风速。
第七步，对飞行器周围风场的测量包含风向角，风向角包含风向方位角和风向高低角，风速与北向的夹角称为风向方位角
Figure BDA0002546877930000034
风速与水平面的夹角称为风向高低角
Figure BDA0002546877930000035
第八步，以惯性导航系统惯性测量单元测得的三轴加速度[ax ay az]T和角速率[pq r]T作为系统的输入量；选取飞行器的机体三轴速度[u v w]T，横滚角φ、俯仰角θ和偏航角ψ和三轴风速[wx wy wz]T作为状态量，即X＝[u v w φ θ ψ wx wy wz]T，并建立无迹卡尔曼滤波器状态方程：Xk＝f(Xk-1,uk-1)+Wk-1，其中Xk和Xk-1分别是k和k-1时刻的状态向量，Wk-1是状态噪声向量，Wk-1由飞行器机体三轴速度变化率
Figure BDA0002546877930000036
三轴姿态角变化率
Figure BDA0002546877930000037
和三轴风速变化率
Figure BDA0002546877930000038
引起；选取卫星定位系统测得的三轴地速[Vx VyVz]T、真空速Va、攻角α和侧滑角β为量测量，即Z＝[Vx Vy Vz Va α β]T，并建立无迹卡尔曼滤波器量测方程Zk＝h(Xk)+Vk，量测噪声Vk由卫星定位系统、大气数据系统的量测噪声引起；根据第一步获得当前时刻即tk时刻的系统输入量，根据第二步和第五步获得当前时刻即tk时刻的量测信息，根据第四步和第五步获得当前时刻即tk时刻的状态量信息，利用无迹卡尔曼滤波器得到tk+1时刻的状态量的最优估计值，从而实现风速风向的在线测量和实时精确估计。
第九步，进行无迹卡尔曼滤波：
8.1)选定滤波初值：
Figure BDA0002546877930000039
Figure BDA00025468779300000310
其中，E指数学期望，X0为状态向量初值，
Figure BDA00025468779300000311
为状态向量初值的均值，P0为状态向量初值的协方差阵
8.2)计算k-1时刻的sigma样本点：
Figure BDA0002546877930000041
Figure BDA0002546877930000042
Figure BDA0002546877930000043
其中，n为状态维数，λ＝α2(n+κ)-n，α∈[10-4,1]，κ＝3-n；
8.3)确定权值：
Figure BDA0002546877930000044
Figure BDA0002546877930000045
Figure BDA0002546877930000046
其中，
Figure BDA0002546877930000047
和Wi (m)为第i个sigma点状态向量的权值，
Figure BDA0002546877930000048
和Wi (c)为第i个sigma点状态向量的协方差阵权值，β为状态分布参数，β＝2；
8.4)计算k时刻的一步预测模型值：
Figure BDA0002546877930000049
Figure BDA00025468779300000410
Figure BDA00025468779300000411
Figure BDA00025468779300000412
为预测状态向量，Pk/k-1为预测状态向量的协方差矩阵，Qk-1为上一个状态向量噪声矩阵；
8.5)计算k时刻的一步预测样本点：
Figure BDA00025468779300000413
Figure BDA00025468779300000414
Figure BDA00025468779300000415
8.6)量测更新：
Figure BDA00025468779300000416
Figure BDA0002546877930000051
Figure BDA0002546877930000052
Figure BDA0002546877930000053
其中，
Figure BDA0002546877930000054
为预测观测向量，
Figure BDA0002546877930000055
为预测向量和观测向量的协方差矩阵，
Figure BDA0002546877930000056
为预测观测量和状态向量的协方差阵，Rk为观测噪声矩阵；
8.7)滤波更新
Figure BDA0002546877930000057
Figure BDA0002546877930000058
Figure BDA0002546877930000059
其中，Kk为增益矩阵，
Figure BDA00025468779300000510
为估计状态向量，Pk为估计状态向量的协方差矩阵。
第十步，无迹卡尔曼滤波器中，状态量的横滚角φ、俯仰角θ和偏航角ψ由惯性导航系统提供，机体系三轴速度由卫星定位系统提供，三维风速初始值有外部风速仪提供，三维风速变化率初始值均取为零，以上组成系统初始状态量，初始状态估计误差方差阵P0，系统噪声方差阵初值Q0，量测噪声方差阵R0矩阵初值都由系统外部直接输入，系统噪声矩阵Qk有系统噪声方差阵初值Q0确定，量测噪声矩阵Rk+1由量测噪声方差阵初值R0确定。
第十一步，由状态量后三项即可得知tk+1时刻的三维风速，将得到的三维风速带到第七步即可得到对应的风向角度。
本发明的优点是：本发明不需要增加额外的机载航电设备，也不需要改变飞行器结构外形，利用现有的惯性导航系统、卫星定位系统和大气数据系统的输出，解算得到飞行器周围的实时风速风向，并能实现飞行航迹上的大区域连续测风，对四维航迹飞行和飞机的准点到达意义重大。
附图说明
图1为一种基于惯性/卫星/大气组合的气象风测量方法原理框图；
图2为一种基于惯性/卫星/大气组合的气象风测量方法的测量方法示意图；
图3为一个实施例中基于惯性/卫星/大气组合的气象风测量方法得到的前向风速估计值和真值误差曲线图；
图4为一个实施例中基于惯性/卫星/大气组合的气象风测量方法得到的右向风速估计值和真值误差曲线图；
图5为一个实施例中基于惯性/卫星/大气组合的气象风测量方法得到的天向风速估计值和真值误差曲线图；
图6为一个实施例中基于惯性/卫星/大气组合的气象风测量方法得到的风向方位角曲线图；
图7为一个实施例中基于惯性/卫星/大气组合的气象风测量方法得到的风向高低角曲线图。
具体实施方式
为了更清楚的说明本申请，下面结合实施例和附图对本申请做进一步的说明。附图中相似的部件以相同的附图标记进行表示。本领域方法人员应当理解，下面所具体描述的内容是说明性的而非限制性的，不应该以此限制本申请的保护范围。
显然，本申请的上述实施例仅仅是为了清楚地说明本发明所作的举例，而非是对本申请的实施方案的限定，对于所属领域的方法人员来说，在上述说明的基础上还可以做出其他不同形式的变化或变动，在这里无法对所有的实施方案予以穷举，凡是属于本申请的方法方案所引申出的显而易见的变化或变动仍处于本发明的保护范围之列。
在一个实施例中，参考图1，本发明利用现有的机载惯性导航系统、卫星定位系统、和大气数据系统的输出参数，通过设计卡尔曼滤波器，实现当前时刻风速风向的在线测量和对下一时刻风速风向的实时精确估计。
为实现风速风向的在线测量和估计，需要以下步骤：
第一步，获取惯性导航系统的输出信息，以周期△T读取惯性导航系统输出的三个角速度、三个姿态角、三个线加速度信息，三个角速度[p q r]T分别是机体系下X轴、Y轴和Z轴方向上的横滚角速度p，俯仰角速度q，偏航角速度r；三个姿态角[φ θ ψ]T分别为横滚角φ、俯仰角θ和偏航角ψ；三个线加速度[ax ay az]T分别为机体坐标系下X轴方向的线加速度ax、Y轴方向的线加速度ay和Z轴方向的线加速度az，其中机体坐标系的X轴、Y轴、Z轴的指向分别为向前、向右和向下。
第二步，获得卫星定位系统的输出信息，以周期ΔT读取机载卫星定位系统输出的地理坐标系下的三轴地速[Vx Vy Vz]T，其中Vx为飞行器在地理坐标系X轴的速度分量、Vy为飞行器在地理坐标系Y轴的速度分量、Vz为飞行器在地理坐标系的速度分量，其中地理坐标系的X轴、Y轴、Z轴的指向分别为北向、东向和地。
第三步，求得系统的状态转移矩阵，状态转移矩阵为地理坐标系到机体坐标系的状态转换矩阵，通过姿态角信息获得，姿态角信息包括横滚角φ、俯仰角θ和偏航角ψ，计算状态转移矩阵
Figure BDA0002546877930000071
Figure BDA0002546877930000072
第四步，通过状态转移矩阵可以地理坐标系下的三轴地速[Vx Vy Vz]T转换到机体系下：
Figure BDA0002546877930000073
[wx wy wz]T为三轴风速。
第五步，获取机载大气数据系统输出数据，机载真空速通常由机载大气数据系统中压力传感器及温度传感器提供的压力以及温度，结合伯努利方程解算得到的，计算公式：
Figure BDA0002546877930000074
其中，△P为动压，ρ为空气密度，K为校正因子。
第六步，利用卫星定位系统测得的地速
Figure BDA0002546877930000075
惯性导航系统测得的姿态角，大气数据系统提供的攻角α和侧滑角β解算得到的转态转移矩阵、大气数据系统测得的真空速
Figure BDA0002546877930000076
得到：
Figure BDA0002546877930000077
其中
Figure BDA0002546877930000081
wx为X轴方向上的风速，wy为Y轴方向上的风速，wz是Z轴方向上的风速。
第七步，从第六步得到三轴风速，代入风向方位角
Figure BDA0002546877930000082
风向高低角
Figure BDA0002546877930000083
求得对应点的风向角。
第八步，以惯性导航系统惯性测量单元测得的三轴加速度[ax ay az]T和角速率[pq r]T作为系统的输入量；选取飞行器的机体三轴速度[u v w]T，横滚角φ、俯仰角θ和偏航角ψ和三轴风速[wx wy wz]T作为状态量，即X＝[u v w φ θ ψ wx wy wz]T，并建立无迹卡尔曼滤波器状态方程：Xk＝f(Xk-1,uk-1)+Wk-1，其中Xk和Xk-1分别是k和k-1时刻的状态向量，Wk-1是状态噪声向量，Wk-1由飞行器机体三轴速度变化率
Figure BDA0002546877930000084
三轴姿态角变化率
Figure BDA0002546877930000085
和三轴风速变化率
Figure BDA0002546877930000086
引起；选取卫星定位系统测得的三轴地速[Vx VyVz]T、真空速Va、攻角α和侧滑角β为量测量，即Z＝[Vx Vy Vz Va α β]T，并建立无迹卡尔曼滤波器量测方程Zk＝h(Xk)+Vk，量测噪声Vk由卫星定位系统、大气数据系统的量测噪声引起；根据第一步获得当前时刻即tk时刻的系统输入量，根据第二步和第五步获得当前时刻即tk时刻的量测信息，根据第四步和第五步获得当前时刻即tk时刻的状态量信息，利用无迹卡尔曼滤波器得到tk+1时刻的状态量的最优估计值，从而实现风速风向的在线测量和实时精确估计。
第九步，进行无迹卡尔曼滤波：
8.1)选定滤波初值：
Figure BDA0002546877930000087
Figure BDA0002546877930000088
其中，E指数学期望，X0为状态向量初值，
Figure BDA0002546877930000089
为状态向量初值的均值，P0为状态向量初值的协方差阵
8.2)计算k-1时刻的sigma样本点：
Figure BDA00025468779300000810
Figure BDA00025468779300000811
Figure BDA0002546877930000091
其中，n为状态维数，λ＝α2(n+κ)-n，α∈[10-4,1]，κ＝3-n；
8.3)确定权值：
Figure BDA0002546877930000092
Figure BDA0002546877930000093
Figure BDA0002546877930000094
其中，
Figure BDA0002546877930000095
和Wi (m)为第i个sigma点状态向量的权值，
Figure BDA0002546877930000096
和Wi (c)为第i个sigma点状态向量的协方差阵权值，β为状态分布参数，β＝2；
8.4)计算k时刻的一步预测模型值：
Figure BDA0002546877930000097
Figure BDA0002546877930000098
Figure BDA0002546877930000099
Figure BDA00025468779300000910
为预测状态向量，Pk/k-1为预测状态向量的协方差矩阵，Qk-1为上一个状态向量噪声矩阵；
8.5)计算k时刻的一步预测样本点：
Figure BDA00025468779300000911
Figure BDA00025468779300000912
Figure BDA00025468779300000913
8.6)量测更新：
Figure BDA00025468779300000914
Figure BDA00025468779300000915
Figure BDA00025468779300000916
Figure BDA0002546877930000101
其中，
Figure BDA0002546877930000102
为预测观测向量，
Figure BDA0002546877930000103
为预测向量和观测向量的协方差矩阵，
Figure BDA0002546877930000104
为预测观测量和状态向量的协方差阵，Rk为观测噪声矩阵；
8.7)滤波更新
Figure BDA0002546877930000105
Figure BDA0002546877930000106
Figure BDA0002546877930000107
其中，Kk为增益矩阵，
Figure BDA0002546877930000108
为估计状态向量，Pk为估计状态向量的协方差矩阵。
第十步，无迹卡尔曼滤波器中，状态量的横滚角φ、俯仰角θ和偏航角ψ由惯性导航系统提供，机体系三轴速度由卫星定位系统提供，三维风速初始值有外部风速仪提供，三维风速变化率初始值均取为零，以上组成系统初始状态量，初始状态估计误差方差阵P0，系统噪声方差阵初值Q0，量测噪声方差阵R0矩阵初值都由系统外部直接输入，系统噪声矩阵Qk有系统噪声方差阵初值Q0确定，量测噪声矩阵Rk+1由量测噪声方差阵初值R0确定。
第十一步，由状态量后三项即可得知tk+1时刻的三维风速，将得到的三维风速带到第七步即可得到tk+1时刻对应的风向角度。
至此，完成该飞行器机体周围的风速风向在线测量和风速风向的实时精确估计。
图1是本发明方法的基本原理框图，无迹卡尔曼滤波器利用机载惯性导航系统、卫星定位系统、大气数据系统的输出，解算得到飞行器机体周围实时风速风向和下一时刻的风速风向估计值。各个模块之间的实线箭头代表了基本的逻辑连接关系。
图2对图1中的原理进行了进一步的细化，代表了基于惯性/卫星/大气组合的气象风测量方法的具体流程，更清晰的表明了系统的数据来源、数据传递和风速风向解算。
图3是本发明方法的一个实施例中得到的前向风速估计值和真值误差曲线图，横轴表示时间，单位为秒，纵轴表示风速估计值与风速真值在前向上的误差，单位为米/秒，范围在[-0.4m/s,0.4m/s]。
图4是本发明方法的一个实施例中得到的右向风速估计值和真值误差曲线图，横轴表示时间，单位为秒，纵轴表示风速估计值与风速真值在右向上的误差，单位为米/秒，范围在[-0.5m/s,0.5m/s]。
图5是本发明方法的一个实施例中得到的天向风速估计值和真值误差曲线图，横轴表示时间，单位为秒，纵轴表示风速估计值与风速真值在天向上的误差，单位为米/秒，范围在[-0.4m/s,0.4m/s]。
图6是本发明方法的一个实施例中得到的风向方位角曲线图，横轴表示时间，单位为秒，纵轴表示风向方位角数值，单位为度，范围在[6°,16°]。
图7是本发明方法的一个实施例中得到的风向高低角曲线图，横轴表示时间，单位为秒，纵轴表示风向高低角数值，单位为度，范围在[17°,27°]。

Claims (10)

1.一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：包括以下步骤：
利用飞行器机载航电设备输出参数，参数包括：
(1)惯性导航系统输出飞行器的加速度、角速度和姿态角；
(2)卫星定位系统输出飞行器的地速；
(3)大气数据系统输出飞行器的真空速；
利用惯性导航系统输出的加速度、角速度和姿态角，对状态转换矩阵进行计算；
利用卫星定位系统输出的地速通过机体坐标系和地理坐标系的状态转移矩阵，求解出机体坐标系下的飞行器三轴速度；
利用大气数据系统输出的真空速，对飞行器的真空速进行解算；
利用当前时刻飞行器的地速、状态转移矩阵和真空速解算结果，对飞行器机身周围的风速风向进行计算，并通过无迹卡尔曼滤波器对下一时刻风速风向进行实时精确估计。
2.根据权利要求1所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的飞行器的加速度、角速度和姿态角的获取方式具体如下：以周期△T读取惯性导航系统输出的三个角速度、三个姿态角和三个线加速度信息，三个角速度[p q r]T分别是机体系下X轴、Y轴和Z轴方向上的横滚角速度p、俯仰角速度q和偏航角速度r；三个姿态角[φ θ ψ]T分别为横滚角φ、俯仰角θ和偏航角ψ；三个线加速度[ax ay az]T分别为机体坐标系下X轴方向的线加速度ax、Y轴方向的线加速度ay和Z轴方向的线加速度az，其中机体坐标系的X轴、Y轴、Z轴的指向分别为向前、向右和向下。
3.根据权利要求2所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的飞行器的低速的获取方式具体如下：以周期△T读取机载卫星定位系统输出的地理坐标系下的三轴地速[Vx Vy Vz]T，其中Vx为飞行器在地理坐标系X轴的速度分量，Vy为飞行器在地理坐标系Y轴的速度分量，Vz为飞行器在地理坐标系的速度分量，其中地理坐标系的X轴、Y轴、Z轴的指向分别为北向、东向和地。
4.根据权利要求3所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的状态转移矩阵为地理坐标系到机体坐标系的状态转换矩阵，通过姿态角信息获得，姿态角信息包括横滚角φ、俯仰角θ和偏航角ψ，计算状态转移矩阵
Figure FDA0002546877920000021
Figure FDA0002546877920000022
通过状态转移矩阵将地理坐标系下的三轴地速[Vx Vy Vz]T转换到机体系下：
Figure FDA0002546877920000023
[wx wy wz]T为三轴风速。
5.根据权利要求4所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的飞行器的真空速由机载大气数据系统中压力传感器及温度传感器提供的压力以及温度，结合伯努利方程解算得到的，计算公式：
Figure FDA0002546877920000024
其中，△P为动压，ρ为空气密度，K为校正因子；
利用卫星定位系统测得的地速
Figure FDA0002546877920000025
惯性导航系统测得的姿态角，大气数据系统提供的攻角α和侧滑角β解算得到的状态转移矩阵、大气数据系统测得的真空速
Figure FDA0002546877920000026
得到
Figure FDA0002546877920000027
其中
Figure FDA0002546877920000028
wx为X轴方向上的风速，wy为Y轴方向上的风速，wz是Z轴方向上的风速。
6.根据权利要求5所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的飞行器机身周围的风向包含风向角，风向角包含风向方位角和风向高低角，风速与北向的夹角称为风向方位角
Figure FDA0002546877920000031
风速与水平面的夹角称为风向高低角
Figure FDA0002546877920000032
7.根据权利要求6所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的对下一时刻风速风向进行实时精确估计，具体方法如下：以三轴加速度[ax ayaz]T和角速率[p q r]T作为系统的输入量；以飞行器的机体三轴速度[u v w]T、横滚角φ、俯仰角θ、偏航角ψ和三轴风速[wx wy wz]T作为状态量，即X＝[u v w φ θ ψ wx wy wz]T；以三轴地速[Vx Vy Vz]T、真空速Va、攻角α和侧滑角β为量测量，即Z＝[Vx Vy Vz Va α β]T；根据获得当前时刻即tk时刻的系统输入量、当前时刻即tk时刻的量测信息和当前时刻即tk时刻的状态量信息，利用无迹卡尔曼滤波器得到tk+1时刻的状态量的最优估计值，从而实现风速风向的在线测量和实时精确估计。
8.根据权利要求7所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的利用无迹卡尔曼滤波器得到tk+1时刻的状态量的最优估计值，从而实现风速风向的在线测量和实时精确估计，具体方法如下：建立无迹卡尔曼滤波器状态方程：Xk＝f(Xk-1,uk-1)+Wk-1，其中Xk和Xk-1分别是k和k-1时刻的状态向量，Wk-1是状态噪声向量，Wk-1由飞行器机体三轴速度变化率
Figure FDA0002546877920000033
三轴姿态角变化率
Figure FDA0002546877920000034
和三轴风速变化率
Figure FDA0002546877920000035
引起；建立无迹卡尔曼滤波器量测方程：Zk＝h(Xk)+Vk，其中Vk为观测噪声向量，Z＝[Vx Vy Vz Va α β]T，量测噪声Vk由卫星定位系统、大气数据系统的量测噪声引起，无迹卡尔曼滤波具体表示为：
8.1)选定滤波初值：
Figure FDA0002546877920000036
Figure FDA0002546877920000037
其中，E指数学期望，X0为状态向量初值，
Figure FDA0002546877920000038
为状态向量初值的均值，P0为状态向量初值的协方差阵；
8.2)计算k-1时刻的sigma样本点：
Figure FDA0002546877920000041
Figure FDA0002546877920000042
Figure FDA0002546877920000043
其中，n为状态维数，λ＝α2(n+κ)-n，α∈[10-4,1]，κ＝3-n；
8.3)确定权值：
Figure FDA0002546877920000044
Figure FDA0002546877920000045
Figure FDA0002546877920000046
其中，W0 (m)和Wi (m)为第i个sigma点状态向量的权值，W0 (c)和Wi (c)为第i个sigma点状态向量的协方差阵权值，β为状态分布参数，β＝2；
8.4)计算k时刻的一步预测模型值：
Figure FDA0002546877920000047
Figure FDA0002546877920000048
Figure FDA0002546877920000049
Figure FDA00025468779200000410
为预测状态向量，Pk/k-1为预测状态向量的协方差矩阵，Qk-1为上一个状态向量噪声矩阵；
8.5)计算k时刻的一步预测样本点
Figure FDA00025468779200000411
Figure FDA00025468779200000412
Figure FDA00025468779200000413
8.6)量测更新：
Figure FDA00025468779200000414
Figure FDA0002546877920000051
Figure FDA0002546877920000052
Figure FDA0002546877920000053
其中，
Figure FDA0002546877920000054
为预测观测向量，
Figure FDA0002546877920000059
为预测向量和观测向量的协方差矩阵，
Figure FDA00025468779200000510
为预测观测量和状态向量的协方差阵，Rk为观测噪声矩阵；
8.7)滤波更新：
Figure FDA0002546877920000055
Figure FDA0002546877920000056
Figure FDA0002546877920000057
其中，Kk为增益矩阵，
Figure FDA0002546877920000058
为估计状态向量，Pk为估计状态向量的协方差矩阵。
9.根据权利要求8所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：所述的无迹卡尔曼滤波器中，状态量的横滚角φ、俯仰角θ和偏航角ψ由惯性导航系统提供，机体系三轴速度由卫星定位系统提供，三维风速初始值由外部风速仪提供，三维风速变化率初始值均取为零，以上组成系统初始状态量，初始状态估计误差方差阵P0、系统噪声方差阵初值Q0和量测噪声方差阵R0矩阵初值都由系统外部直接输入，系统噪声矩阵Qk由系统噪声方差阵初值Q0确定，量测噪声矩阵Rk+1由量测噪声方差阵初值R0确定。
10.根据权利要求9所述的一种基于惯性/卫星/大气组合的气象风测量方法，其特征在于：由所述的状态量得知tk+1时刻的三维风速，将得到的三维风速带到风速与北向的夹角公式和风速与水平面的夹角公式中得到风向角度。
CN202010563456.4A 2020-06-19 2020-06-19 一种基于惯性/卫星/大气组合的气象风测量方法 Active CN111766397B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010563456.4A CN111766397B (zh) 2020-06-19 2020-06-19 一种基于惯性/卫星/大气组合的气象风测量方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010563456.4A CN111766397B (zh) 2020-06-19 2020-06-19 一种基于惯性/卫星/大气组合的气象风测量方法

Publications (2)

Publication Number Publication Date
CN111766397A true CN111766397A (zh) 2020-10-13
CN111766397B CN111766397B (zh) 2022-06-10

Family

ID=72721456

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010563456.4A Active CN111766397B (zh) 2020-06-19 2020-06-19 一种基于惯性/卫星/大气组合的气象风测量方法

Country Status (1)

Country Link
CN (1) CN111766397B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111965667A (zh) * 2020-10-14 2020-11-20 南京牧镭激光科技有限公司 动态补偿测风激光雷达系统及其测风方法
CN114636842A (zh) * 2022-05-17 2022-06-17 成都信息工程大学 一种高超声速飞行器的大气数据估计方法及装置
CN117251942A (zh) * 2023-11-17 2023-12-19 成都凯天电子股份有限公司 一种估算飞行器空速、攻角和侧滑角的方法及系统

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20040024500A1 (en) * 2002-07-30 2004-02-05 Campbell Charles R. Method for estimating wind
CN102565451A (zh) * 2011-12-28 2012-07-11 中国航空工业集团公司第六三一研究所 通用飞机航行风向、风速的测算方法
CN103852081A (zh) * 2014-03-20 2014-06-11 南京航空航天大学 用于大气数据/捷联惯导组合导航系统的真空速解算方法
CN104459193A (zh) * 2014-12-05 2015-03-25 中国航天空气动力技术研究院 一种基于无人机侧航法估算侧风信息的方法
CN105388535A (zh) * 2015-11-11 2016-03-09 上海埃威航空电子有限公司 基于现有机载设备的航空气象风观测方法
CN106885918A (zh) * 2017-02-10 2017-06-23 南京航空航天大学 一种面向多旋翼飞行器的多信息融合实时风速估计方法
KR101862065B1 (ko) * 2017-07-25 2018-05-29 한국항공대학교산학협력단 비행체를 이용한 영상 기반 바람 추정 장치 및 방법
CN110873813A (zh) * 2019-12-02 2020-03-10 中国人民解放军战略支援部队信息工程大学 一种水流速度估算方法、组合导航方法及装置

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20040024500A1 (en) * 2002-07-30 2004-02-05 Campbell Charles R. Method for estimating wind
CN102565451A (zh) * 2011-12-28 2012-07-11 中国航空工业集团公司第六三一研究所 通用飞机航行风向、风速的测算方法
CN103852081A (zh) * 2014-03-20 2014-06-11 南京航空航天大学 用于大气数据/捷联惯导组合导航系统的真空速解算方法
CN104459193A (zh) * 2014-12-05 2015-03-25 中国航天空气动力技术研究院 一种基于无人机侧航法估算侧风信息的方法
CN105388535A (zh) * 2015-11-11 2016-03-09 上海埃威航空电子有限公司 基于现有机载设备的航空气象风观测方法
CN106885918A (zh) * 2017-02-10 2017-06-23 南京航空航天大学 一种面向多旋翼飞行器的多信息融合实时风速估计方法
KR101862065B1 (ko) * 2017-07-25 2018-05-29 한국항공대학교산학협력단 비행체를 이용한 영상 기반 바람 추정 장치 및 방법
CN110873813A (zh) * 2019-12-02 2020-03-10 中国人民解放军战略支援部队信息工程大学 一种水流速度估算方法、组合导航方法及装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
何波等: "《基于扩展卡尔曼滤波的风速估计算法研究》", 《电子测量技术》 *
曹良秋等: "《基于无人机飞行参数的风场估算技术研究》", 《科技创新导报》 *
李静等: "《基于气动角的气象无人机风场测量和数据处理方法》", 《智能计算机与应用》 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111965667A (zh) * 2020-10-14 2020-11-20 南京牧镭激光科技有限公司 动态补偿测风激光雷达系统及其测风方法
CN111965667B (zh) * 2020-10-14 2020-12-29 南京牧镭激光科技有限公司 动态补偿测风激光雷达系统及其测风方法
US11821991B2 (en) 2020-10-14 2023-11-21 Nanjing Movelaser Co., Ltd. Dynamic compensation wind measurement lidar system and wind measurement method thereof
CN114636842A (zh) * 2022-05-17 2022-06-17 成都信息工程大学 一种高超声速飞行器的大气数据估计方法及装置
CN117251942A (zh) * 2023-11-17 2023-12-19 成都凯天电子股份有限公司 一种估算飞行器空速、攻角和侧滑角的方法及系统
CN117251942B (zh) * 2023-11-17 2024-03-08 成都凯天电子股份有限公司 一种估算飞行器空速、攻角和侧滑角的方法及系统

Also Published As

Publication number Publication date
CN111766397B (zh) 2022-06-10

Similar Documents

Publication Publication Date Title
CN111766397B (zh) 一种基于惯性/卫星/大气组合的气象风测量方法
CN109813311B (zh) 一种无人机编队协同导航方法
CN111024064B (zh) 一种改进Sage-Husa自适应滤波的SINS/DVL组合导航方法
CN109033485B (zh) 用于基于天气缓冲模型估计飞行器空速的系统
CN107270893B (zh) 面向不动产测量的杆臂、时间不同步误差估计与补偿方法
CN110780326A (zh) 一种车载组合导航系统和定位方法
CN113124856B (zh) 基于uwb在线锚点的视觉惯性紧耦合里程计及计量方法
CN110567454B (zh) 一种复杂环境下sins/dvl紧组合导航方法
CN112505737B (zh) 一种gnss/ins组合导航方法
CN101122637A (zh) 一种sar运动补偿用sins/gps组合导航自适应降维滤波方法
CN112556724A (zh) 动态环境下的微型飞行器低成本导航系统初始粗对准方法
Youn et al. Aerodynamic model-aided estimation of attitude, 3-D wind, airspeed, AOA, and SSA for high-altitude long-endurance UAV
CN110007318B (zh) 风场干扰下基于卡尔曼滤波的单无人机判断gps欺骗的方法
CN111964688A (zh) 结合无人机动力学模型和mems传感器的姿态估计方法
CN116992700B (zh) 一种物流无人机导航精度确定的方法及设备
CN114636842B (zh) 一种高超声速飞行器的大气数据估计方法及装置
Lu et al. Air data estimation by fusing navigation system and flight control system
CN113008229A (zh) 一种基于低成本车载传感器的分布式自主组合导航方法
CN115542363B (zh) 一种适用于垂直下视航空吊舱的姿态测量方法
CN116576849A (zh) 一种基于gmm辅助的车辆融合定位方法及系统
Emran et al. A cascaded approach for quadrotor's attitude estimation
Yuan et al. Dynamic initial alignment of the MEMS-based low-cost SINS for AUV based on unscented Kalman filter
CN111308127B (zh) 一种基于大气物理机制的星载加速度计校准方法
CN112762960A (zh) 一种飞行器所处风场的在线计算方法
Gong et al. Unscented particle smoother and its application to transfer alignment of airborne distributed POS

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant