WO2020220729A1 - 基于角加速度计/陀螺/加速度计的惯性导航解算方法 - Google Patents

基于角加速度计/陀螺/加速度计的惯性导航解算方法 Download PDF

Info

Publication number
WO2020220729A1
WO2020220729A1 PCT/CN2020/070025 CN2020070025W WO2020220729A1 WO 2020220729 A1 WO2020220729 A1 WO 2020220729A1 CN 2020070025 W CN2020070025 W CN 2020070025W WO 2020220729 A1 WO2020220729 A1 WO 2020220729A1
Authority
WO
WIPO (PCT)
Prior art keywords
time
accelerometer
carrier
angular
angular velocity
Prior art date
Application number
PCT/CN2020/070025
Other languages
English (en)
French (fr)
Inventor
吕品
白师宇
赖际舟
李志敏
刘士超
Original Assignee
南京航空航天大学
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 南京航空航天大学 filed Critical 南京航空航天大学
Publication of WO2020220729A1 publication Critical patent/WO2020220729A1/zh

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G01C21/10Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
    • G01C21/12Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
    • G01C21/16Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
    • G01C21/18Stabilised platforms, e.g. by gyroscope
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C25/00Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C25/00Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass
    • G01C25/005Manufacturing, calibrating, cleaning, or repairing instruments or devices referred to in the other groups of this subclass initial alignment, calibration or starting-up of inertial devices
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Definitions

  • the invention relates to an inertial navigation solution method based on an angular accelerometer/gyro/accelerometer, and belongs to the technical field of inertial navigation.
  • Inertial navigation is a common navigation method that uses inertial devices to calculate the attitude, speed, and position of the carrier through a recursive navigation method. Inertial navigation has the advantages of strong autonomy, freedom from external interference, and complete output information. It has a wide range of applications in aviation, aerospace, and navigation.
  • the inertial devices used in traditional inertial navigation solutions are gyroscopes and accelerometers.
  • the gyroscope measures angular velocity information
  • the accelerometer measures acceleration information.
  • the bandwidth of the gyro is limited, and the sampling delay will bring errors, thereby affecting the navigation accuracy.
  • the technical problem to be solved by the present invention is to provide an inertial navigation solution method based on angular accelerometer/gyro/accelerometer, form a new inertial navigation scheme by introducing angular accelerometer, and design corresponding inertial navigation algorithm to improve The navigation accuracy of the carrier in a highly dynamic environment.
  • the inertial navigation solution method based on angular accelerometer/gyro/accelerometer includes the following steps:
  • Step 1 Periodically read the sensor information of the carrier at time k, including the output of the three-axis angular accelerometer Three-axis gyro output Three-axis accelerometer output
  • Step 2 Predict the angular velocity, attitude, velocity, and position information of the carrier at time k based on the sensor information at time k;
  • Step 3 Correct the angular velocity of the carrier at time k through the Kalman filter.
  • step 2 the specific process of step 2 is as follows:
  • q 0 (k-1), q 1 (k-1), q 2 (k-1), q 3 (k-1) are the attitude quaternions at time k-1
  • Is the attitude transition matrix at time k-1 Is the component of the angular velocity of the navigation system relative to the inertial system on the X, Y, and Z axes of the navigation system at k-1
  • Is the attitude transition matrix at time k-1 Is the component of the angular velocity of the navigation system relative to the inertial system on the X, Y, and Z axes of the navigation system at k-1
  • Is the component of the carrier velocity on the X, Y, and Z axes of the navigation system at time k Is the component of the carrier speed on the X, Y, and Z axes of the navigation system at time k-1, L(k-1), h(k-1) are the latitude and height of the carrier at time k-1, and R M and R N are The radius of the earth’s meridian circle and Maoyou circle, Is the output component of the three-axis accelerometer on the X, Y, and Z axes at time k, ⁇ ie is the angular velocity of the earth's rotation, and g is the acceleration of gravity;
  • ⁇ (k), L(k), h(k) are the longitude, latitude, and height of the carrier at time k
  • ⁇ (k-1) is the longitude of the carrier at time k-1.
  • the attitude transition matrix at time k-1 for:
  • q 0 (k-1), q 1 (k-1), q 2 (k-1), and q 3 (k-1) are posture quaternions at time k-1.
  • the angular velocity of the navigation system relative to the inertial system at the time k-1 is the component of the navigation system on the X, Y, and Z axes for:
  • step 3 the specific process of step 3 is as follows:
  • k-1) A(k,k-1)P(k-1
  • A(k,k-1) is the filter step transition matrix from time k-1 to time k
  • k-1) is the mean square error of state estimation at time k-1
  • k-1) is the one-step prediction mean square error from time k-1 to time k
  • G(k-1) is the filter noise coefficient matrix at time k-1
  • W(k-1) is k-1 State noise at the moment
  • T means transpose
  • K(k) P(k
  • K(k) is the filter gain at time k
  • H(k) is the measurement matrix at time k
  • R(k) is the measurement noise at time k
  • superscript -1 means inversion
  • Is the estimated value of the state quantity at time k Is the one-step predicted value of the state variable from k-1 to k, Y(k) is the measured value at k, Is the output component of the three-axis gyroscope on the X, Y, and Z axes at time k;
  • k) is the estimated mean square error at time k
  • I is the identity matrix
  • the present invention adopts the above technical solutions and has the following technical effects:
  • the invention can reduce the angular velocity estimation error caused by the limited bandwidth of the gyro and the sampling delay of the carrier in the high dynamic environment, and improve the angular velocity and attitude estimation accuracy.
  • Fig. 1 is a schematic flowchart of the inertial navigation solution method based on angular accelerometer/gyro/accelerometer of the present invention.
  • Fig. 2, Fig. 3, and Fig. 4 are respectively the X-axis, Y-axis and Z-axis angular velocity errors of the aircraft in the flight state when the traditional inertial calculation method and the method of the present invention are used.
  • Fig. 5, Fig. 6, Fig. 7 are respectively the roll angle, pitch angle and heading angle error of the aircraft when the traditional inertial calculation method and the method of the present invention are used.
  • Step 1 Periodically read sensor information at time k, including three-axis angular accelerometer Three-axis gyro Three-axis accelerometer
  • Step 2 Predict the angular velocity, attitude, velocity, and position information of the carrier at time k:
  • q 0 (k), q 1 (k), q 2 (k), q 3 (k) are the attitude quaternions at time k
  • q 0 (k-1), q 1 (k-1), q 2 (k-1), q 3 (k-1) are the posture quaternions at time k-1
  • Is the attitude transition matrix at time k-1 Is the component of the angular velocity of the navigation system relative to the inertial system on the X, Y, and Z axes of the navigation system at k-1
  • L(k-1), h(k-1) are the latitude and height of the carrier at k-1
  • R M and R N are The radius of the earth’s meridian and Maoyou
  • ⁇ (k), L(k), h(k) are the longitude, latitude and height at time k
  • ⁇ (k-1), L(k-1), h(k-1) are time k-1 Longitude, latitude and altitude.
  • Step 3 Use Kalman filter to correct the angular velocity of the carrier at time k:
  • k-1) A(k,k-1)P(k-1
  • A(k,k-1) is the one-step transition matrix of the filter from time k-1 to time k
  • I 3 ⁇ 3 is a 3 ⁇ 3 unit matrix
  • 0 3 ⁇ 3 is a 3 ⁇ 3 zero matrix
  • 0 3 ⁇ 4 is a 3 ⁇ 4 zero matrix
  • 0 4 ⁇ 3 is a 4 ⁇ 3 zero matrix
  • k-1) is the mean square error of state estimation at time k-1
  • k-1) is the one-step predicted mean square error from time k-1 to time k
  • G(k-1) is The filter noise coefficient matrix at the time of filter k-1
  • T is the state noise at time k-1, ⁇ rx , ⁇ ry and ⁇ rz are respectively with The model noise of ⁇ ax , ⁇ ay and ⁇ az are respectively with The model noise of, the superscript T means transpose;
  • K(k) P(k
  • H(k) [I 3 ⁇ 3 0 3 ⁇ 4 0 3 ⁇ 3 0 3 ⁇ 3 ]
  • H(k) is the measurement matrix at k time
  • K(k) is the filter gain at k time
  • diag represents matrix diagonalization, where ⁇ ⁇ x , ⁇ ⁇ y , ⁇ ⁇ z are respectively Noise, superscript -1 means inversion;
  • I the estimated value of the state quantity at time k
  • k) is the estimated mean square error at time k
  • I is the identity matrix
  • the navigation accuracy of the aircraft after using the method of the present invention is verified.
  • the accuracy of the inertial sensor is set as follows: the zero-bias stability of the angular accelerometer is 150deg/s 2 , the zero-bias stability of the gyro is 0.01deg/h, and the zero-bias stability of the accelerometer is 30 ⁇ g.
  • the sampling frequency of angular accelerometer, gyroscope and accelerometer is 400Hz, and the cutoff frequency of gyroscope is 100Hz.
  • Fig. 2, Fig. 3 and Fig. 4 show the angular velocity error of the aircraft in the trajectory flight state when the traditional method and the method of the present invention are used.
  • Fig. 5, Fig. 6 and Fig. 7 show the attitude error of the aircraft under the flight path when the traditional method and the method of the present invention are used.

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Manufacturing & Machinery (AREA)
  • Automation & Control Theory (AREA)
  • Navigation (AREA)

Abstract

一种基于角加速度计/陀螺/加速度计的惯性导航解算方法,其采用角加速度计、陀螺、加速度计传感器对载体的角速度、姿态、速度、位置信息进行预测,通过卡尔曼滤波器对载体的角速度、姿态、速度、位置信息进行解算。其中,角加速度计、加速度计用于更新状态方程,陀螺用于更新量测方程。相对于传统的惯性导航方案,通过引入角加速度计,可以提高载体高动态情况下的导航解算精度。

Description

基于角加速度计/陀螺/加速度计的惯性导航解算方法 技术领域
本发明涉及基于角加速度计/陀螺/加速度计的惯性导航解算方法,属于惯性导航技术领域。
背景技术
惯性导航是一种常见的导航方式,其采用惯性器件,通过递推式的导航方式,对载体的姿态、速度、位置进行解算。惯性导航具有自主性强、不受外界干扰、输出信息全的优势,在航空、航天、航海中具有广泛的应用。
传统的惯性导航方案采用的惯性器件为陀螺与加速度计,其中陀螺测量的是角速度信息,加速度计测量的是加速度信息。当载体进行高动态动作时,陀螺的带宽受限、采样延迟会带来误差,从而影响导航精度。
发明内容
本发明所要解决的技术问题是:提供基于角加速度计/陀螺/加速度计的惯性导航解算方法,通过引入角加速度计,形成一种新的惯性导航方案,并设计相应的惯性导航算法,提高载体在高动态环境下的导航精度。
本发明为解决上述技术问题采用以下技术方案:
基于角加速度计/陀螺/加速度计的惯性导航解算方法,包括如下步骤:
步骤1,周期读取k时刻载体的传感器信息,包括三轴角加速度计输出
Figure PCTCN2020070025-appb-000001
三轴陀螺输出
Figure PCTCN2020070025-appb-000002
三轴加速度计输出
Figure PCTCN2020070025-appb-000003
步骤2,根据k时刻的传感器信息,预测k时刻载体的角速度、姿态、速度、位置信息;
步骤3,通过卡尔曼滤波器,对k时刻载体的角速度进行校正。
作为本发明的一种优选方案,所述步骤2的具体过程如下:
1)采用如下公式预测载体角速度:
Figure PCTCN2020070025-appb-000004
Figure PCTCN2020070025-appb-000005
Figure PCTCN2020070025-appb-000006
其中,
Figure PCTCN2020070025-appb-000007
为k时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000008
为k-1时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000009
为k时刻机体系相对于惯性系的角加速度在机体系X、Y、Z轴上的分量,ΔT为离散采样周期;
2)采用如下公式预测载体姿态:
Figure PCTCN2020070025-appb-000010
Figure PCTCN2020070025-appb-000011
Figure PCTCN2020070025-appb-000012
Figure PCTCN2020070025-appb-000013
其中,
Figure PCTCN2020070025-appb-000014
a=x,y,z,q 0(k)、q 1(k)、q 2(k)、q 3(k)为k时刻的姿态四元数,q 0(k-1)、q 1(k-1)、q 2(k-1)、q 3(k-1)为k-1时刻的姿态四元数,
Figure PCTCN2020070025-appb-000015
为k时刻机体系相对于导航系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000016
为k-1时刻的姿态转移矩阵,
Figure PCTCN2020070025-appb-000017
为k-1时刻导航系相对于惯性系的角速度在导航系X、Y、Z轴上的分量;
3)采用如下公式预测载体速度:
Figure PCTCN2020070025-appb-000018
Figure PCTCN2020070025-appb-000019
Figure PCTCN2020070025-appb-000020
其中,
Figure PCTCN2020070025-appb-000021
为k时刻载体速度在导航系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000022
为k-1时刻载体速度在导航系X、Y、Z轴上的分量,L(k-1)、h(k-1)为k-1时刻载体的纬度、高度,R M、R N为地球的子午圈、卯酉圈半径,
Figure PCTCN2020070025-appb-000023
为k时刻三轴加速度计在X、Y、Z轴上的输出分量,ω ie为地球自转角速度,g为重力加速度;
4)采用如下公式预测载***置:
Figure PCTCN2020070025-appb-000024
Figure PCTCN2020070025-appb-000025
Figure PCTCN2020070025-appb-000026
其中,λ(k)、L(k)、h(k)为k时刻载体的经度、纬度与高度,λ(k-1)为k-1时刻载体的 经度。
作为本发明的一种优选方案,所述k-1时刻的姿态转移矩阵
Figure PCTCN2020070025-appb-000027
为:
Figure PCTCN2020070025-appb-000028
其中,q 0(k-1)、q 1(k-1)、q 2(k-1)、q 3(k-1)为k-1时刻的姿态四元数。
作为本发明的一种优选方案,所述k-1时刻导航系相对于惯性系的角速度在导航系X、Y、Z轴上的分量
Figure PCTCN2020070025-appb-000029
为:
Figure PCTCN2020070025-appb-000030
其中,
Figure PCTCN2020070025-appb-000031
为k-1时刻载体速度在导航系X、Y轴上的分量,L(k-1)、h(k-1)为k-1时刻载体的纬度、高度,R M、R N为地球的子午圈、卯酉圈半径,ω ie为地球自转角速度,上标T表示转置。
作为本发明的一种优选方案,所述步骤3的具体过程如下:
1)计算一步预测均方误差P(k|k-1):
P(k|k-1)=A(k,k-1)P(k-1|k-1)A(k,k-1) T+G(k-1)W(k-1)G(k-1) T
其中,A(k,k-1)为滤波器k-1时刻到k时刻的滤波器一步转移矩阵,P(k-1|k-1)为k-1时刻的状态估计均方差,P(k|k-1)为k-1时刻到k时刻的一步预测均方差,G(k-1)为滤波器k-1时刻的滤波器噪声系数矩阵,W(k-1)为k-1时刻状态噪声,上标T表示转置;
2)计算k时刻扩展卡尔曼滤波器滤波增益K(k):
K(k)=P(k|k-1)H(k) T[H(k)P(k|k-1)H(k) T+R(k)] -1
其中,K(k)为k时刻的滤波增益,H(k)为k时刻量测矩阵,R(k)为k时刻的量测噪声,上标-1表示求逆;
3)计算k时刻扩展卡尔曼滤波器状态估计值
Figure PCTCN2020070025-appb-000032
Figure PCTCN2020070025-appb-000033
其中,
Figure PCTCN2020070025-appb-000034
为k时刻状态量的估计值,
Figure PCTCN2020070025-appb-000035
为k-1到k时刻的状态变量一步预测值,Y(k)为k时刻的量测值,
Figure PCTCN2020070025-appb-000036
Figure PCTCN2020070025-appb-000037
为k时刻三轴陀螺在X、Y、Z轴上的输出分量;
4)计算k时刻扩展卡尔曼滤波器估计均方误差P(k|k):
P(k|k)=[I-K(k)H(k)]P(k|k-1)
其中,P(k|k)为k时刻估计均方误差,I为单位矩阵;
5)基于卡尔曼滤波器,通过三轴陀螺的输出Y(k)对状态中的角速度
Figure PCTCN2020070025-appb-000038
进行校正。
本发明采用以上技术方案与现有技术相比,具有以下技术效果:
本发明可以降低载体在高动态环境下由于陀螺带宽受限、采样延迟引起的角速度 估计误差,提高角速度、姿态估计精度。
附图说明
图1是本发明基于角加速度计/陀螺/加速度计的惯性导航解算方法的流程示意图。
图2、图3、图4分别是采用传统惯性解算方法和本发明方法时飞行器在航迹飞行状态下的X轴、Y轴、Z轴角速度误差。
图5、图6、图7分别是采用传统惯性解算方法和本发明方法时飞行器在航迹飞行状态下的横滚角、俯仰角、航向角误差。
具体实施方式
下面详细描述本发明的实施方式,所述实施方式的示例在附图中示出。下面通过参考附图描述的实施方式是示例性的,仅用于解释本发明,而不能解释为对本发明的限制。
如图1所示,为本发明方法的流程图,其具体步骤如下:
步骤一:周期读取k时刻的传感器信息,包括三轴角加速度计
Figure PCTCN2020070025-appb-000039
Figure PCTCN2020070025-appb-000040
三轴陀螺
Figure PCTCN2020070025-appb-000041
三轴加速度计
Figure PCTCN2020070025-appb-000042
步骤二:预测k时刻载体的角速度、姿态、速度、位置信息:
1)角速度预测采用如下公式:
Figure PCTCN2020070025-appb-000043
Figure PCTCN2020070025-appb-000044
Figure PCTCN2020070025-appb-000045
其中,
Figure PCTCN2020070025-appb-000046
为k时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000047
为k-1时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000048
为k时刻机体系相对于惯性系的角加速度在机体系X、Y、Z轴上的分量,ΔT为离散采样周期;
2)姿态预测采用如下公式:
Figure PCTCN2020070025-appb-000049
Figure PCTCN2020070025-appb-000050
Figure PCTCN2020070025-appb-000051
Figure PCTCN2020070025-appb-000052
式中,
Figure PCTCN2020070025-appb-000053
Figure PCTCN2020070025-appb-000054
Figure PCTCN2020070025-appb-000055
q 0(k)、q 1(k)、q 2(k)、q 3(k)为k时刻的姿态四元数,q 0(k-1)、q 1(k-1)、q 2(k-1)、q 3(k-1)为k-1时刻的姿态四元数,
Figure PCTCN2020070025-appb-000056
为k时刻机体系相对于导航系的角速度在机体系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000057
为k-1时刻的姿态转移矩阵,
Figure PCTCN2020070025-appb-000058
为k-1时刻导航系相对于惯性系的角速度在导航系X、Y、Z轴上的分量,
Figure PCTCN2020070025-appb-000059
为k-1时刻载体速度在导航系X、Y、Z轴上的分量,L(k-1)、h(k-1)为k-1时刻载体的纬度和高度,R M、R N为地球的子午圈与卯酉圈半径,ω ie为地球自转角速度;
3)速度预测采用如下公式:
Figure PCTCN2020070025-appb-000060
Figure PCTCN2020070025-appb-000061
Figure PCTCN2020070025-appb-000062
其中,
Figure PCTCN2020070025-appb-000063
为k时刻载体速度在导航系X、Y、Z轴上的分量;
4)位置预测采用如下公式:
Figure PCTCN2020070025-appb-000064
Figure PCTCN2020070025-appb-000065
Figure PCTCN2020070025-appb-000066
其中,λ(k)、L(k)、h(k)为k时刻经度、纬度与高度,λ(k-1)、L(k-1)、h(k-1)为k-1时刻经度、纬度与高度。
步骤三:通过卡尔曼滤波器,对k时刻载体的角速度进行校正:
1)计算一步预测均方误差P(k|k-1):
P(k|k-1)=A(k,k-1)P(k-1|k-1)A(k,k-1) T+G(k-1)W(k-1)G(k-1) T
式中,A(k,k-1)为滤波器k-1时刻到k时刻的滤波器一步转移矩阵,
Figure PCTCN2020070025-appb-000067
I 3×3为3×3的单位矩阵,0 3×3为3×3的零矩阵,0 3×4为3×4的零矩阵,0 4×3为4×3的零矩阵,
Figure PCTCN2020070025-appb-000068
Figure PCTCN2020070025-appb-000069
Figure PCTCN2020070025-appb-000070
Figure PCTCN2020070025-appb-000071
Figure PCTCN2020070025-appb-000072
P(k-1|k-1)为k-1时刻的状态估计均方差,P(k|k-1)为k-1时刻到k时刻的一步预测均方差,G(k-1)为滤波器k-1时刻的滤波器噪声系数矩阵,
Figure PCTCN2020070025-appb-000073
Figure PCTCN2020070025-appb-000074
Figure PCTCN2020070025-appb-000075
W=[ε rx ε ry ε rz ε ax ε ay ε az] T为k-1时刻状态噪声,ε rx、ε ry和ε rz分别为
Figure PCTCN2020070025-appb-000076
Figure PCTCN2020070025-appb-000077
的模型噪声,ε ax、ε ay和ε az分别为
Figure PCTCN2020070025-appb-000078
Figure PCTCN2020070025-appb-000079
的模型噪声,上标T表示转置;
2)计算k时刻扩展卡尔曼滤波器滤波增益K(k):
K(k)=P(k|k-1)H(k) T[H(k)P(k|k-1)H(k) T+R(k)] -1
式中,H(k)=[I 3×3 0 3×4 0 3×3 0 3×3],H(k)为k时刻量测矩阵,K(k)为k时刻的滤波增益,R(k)=diag([ε ωx ε ωy ε ωz] 2)为k时刻的量测噪声,diag表示矩阵对角化,其中ε ωx、ε ωy、ε ωz分别为
Figure PCTCN2020070025-appb-000080
的噪声,上标-1表示求逆;
3)计算k时刻扩展卡尔曼滤波器状态估计值
Figure PCTCN2020070025-appb-000081
Figure PCTCN2020070025-appb-000082
式中,
Figure PCTCN2020070025-appb-000083
为k时刻状态量的估计值,
Figure PCTCN2020070025-appb-000084
为k-1到k时刻的状态变量一步预测值,使用步骤二的预测公式计算得到,
Figure PCTCN2020070025-appb-000085
为k时刻的量测值,通过步骤一的传感器数据读取获得;
4)计算k时刻扩展卡尔曼滤波器估计均方误差P(k|k):
P(k|k)=[I-K(k)H(k)]P(k|k-1)
式中,P(k|k)为k时刻估计均方误差,I为单位矩阵。
5)基于卡尔曼滤波器,通过陀螺的输出Y(k)对状态中的角速度
Figure PCTCN2020070025-appb-000086
进行校正。
实施例:
采用仿真的形式,对使用本发明方法后的飞行器导航精度进行验证。其中惯性传感器精度设置如下:角加速度计零偏稳定性为150deg/s 2,陀螺零偏稳定性为0.01deg/h,加速度计零偏稳定性为30μg。角加速度计、陀螺与加速度计采样频率均为400Hz,陀螺截止频率为100Hz。
令飞行器从地面加速滑跑起飞,在空中依次分别绕飞机的三个轴做大机动旋转运动,依次为各0.5s的正负90°横滚运动,各0.5s的正负45°俯仰运动,0.5s的正向90°航向运动。
图2、图3以及图4为采用传统方法、本发明方法时飞行器在航迹飞行状态下的角速度误差。
图5、图6以及图7为采用传统方法、本发明方法时飞行器在航迹飞行状态下的姿态误差。
以上实施例仅为说明本发明的技术思想,不能以此限定本发明的保护范围,凡是按照本发明提出的技术思想,在技术方案基础上所做的任何改动,均落入本发明保护范围之内。

Claims (5)

  1. 基于角加速度计/陀螺/加速度计的惯性导航解算方法,其特征在于,包括如下步骤:
    步骤1,周期读取k时刻载体的传感器信息,包括三轴角加速度计输出
    Figure PCTCN2020070025-appb-100001
    三轴陀螺输出
    Figure PCTCN2020070025-appb-100002
    三轴加速度计输出
    Figure PCTCN2020070025-appb-100003
    步骤2,根据k时刻的传感器信息,预测k时刻载体的角速度、姿态、速度、位置信息;
    步骤3,通过卡尔曼滤波器,对k时刻载体的角速度进行校正。
  2. 根据权利要求1所述基于角加速度计/陀螺/加速度计的惯性导航解算方法,其特征在于,所述步骤2的具体过程如下:
    1)采用如下公式预测载体角速度:
    Figure PCTCN2020070025-appb-100004
    Figure PCTCN2020070025-appb-100005
    Figure PCTCN2020070025-appb-100006
    其中,
    Figure PCTCN2020070025-appb-100007
    为k时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
    Figure PCTCN2020070025-appb-100008
    为k-1时刻机体系相对于惯性系的角速度在机体系X、Y、Z轴上的分量,
    Figure PCTCN2020070025-appb-100009
    为k时刻机体系相对于惯性系的角加速度在机体系X、Y、Z轴上的分量,ΔT为离散采样周期;
    2)采用如下公式预测载体姿态:
    Figure PCTCN2020070025-appb-100010
    Figure PCTCN2020070025-appb-100011
    Figure PCTCN2020070025-appb-100012
    Figure PCTCN2020070025-appb-100013
    其中,
    Figure PCTCN2020070025-appb-100014
    a=x,y,z,q 0(k)、q 1(k)、q 2(k)、q 3(k)为k时刻的姿态四元数,q 0(k-1)、q 1(k-1)、q 2(k-1)、q 3(k-1)为k-1时刻的姿态四元数,
    Figure PCTCN2020070025-appb-100015
    为k时刻机体系相对于导航系的角速度在机体系X、Y、Z轴上的分量,
    Figure PCTCN2020070025-appb-100016
    为k-1时刻的姿态转移矩阵,
    Figure PCTCN2020070025-appb-100017
    为k-1时刻导航系相对于惯性系的角速度在导航系X、Y、Z轴上的分量;
    3)采用如下公式预测载体速度:
    Figure PCTCN2020070025-appb-100018
    Figure PCTCN2020070025-appb-100019
    Figure PCTCN2020070025-appb-100020
    其中,
    Figure PCTCN2020070025-appb-100021
    为k时刻载体速度在导航系X、Y、Z轴上的分量,
    Figure PCTCN2020070025-appb-100022
    为k-1时刻载体速度在导航系X、Y、Z轴上的分量,L(k-1)、h(k-1)为k-1时刻载体的纬度、高度,R M、R N为地球的子午圈、卯酉圈半径,
    Figure PCTCN2020070025-appb-100023
    为k时刻三轴加速度计在X、Y、Z轴上的输出分量,ω ie为地球自转角速度,g为重力加速度;
    4)采用如下公式预测载***置:
    Figure PCTCN2020070025-appb-100024
    Figure PCTCN2020070025-appb-100025
    Figure PCTCN2020070025-appb-100026
    其中,λ(k)、L(k)、h(k)为k时刻载体的经度、纬度与高度,λ(k-1)为k-1时刻载体的经度。
  3. 根据权利要求2所述基于角加速度计/陀螺/加速度计的惯性导航解算方法,其特征在于,所述k-1时刻的姿态转移矩阵
    Figure PCTCN2020070025-appb-100027
    为:
    Figure PCTCN2020070025-appb-100028
    其中,q 0(k-1)、q 1(k-1)、q 2(k-1)、q 3(k-1)为k-1时刻的姿态四元数。
  4. 根据权利要求2所述基于角加速度计/陀螺/加速度计的惯性导航解算方法,其特征在于,所述k-1时刻导航系相对于惯性系的角速度在导航系X、Y、Z轴上的分量
    Figure PCTCN2020070025-appb-100029
    为:
    Figure PCTCN2020070025-appb-100030
    其中,
    Figure PCTCN2020070025-appb-100031
    为k-1时刻载体速度在导航系X、Y轴上的分量,L(k-1)、h(k-1)为k-1时刻载体的纬度、高度,R M、R N为地球的子午圈、卯酉圈半径,ω ie为地球自转角速度,上标T表示转置。
  5. 根据权利要求1所述基于角加速度计/陀螺/加速度计的惯性导航解算方法,其特征在于,所述步骤3的具体过程如下:
    1)计算一步预测均方误差P(k|k-1):
    P(k|k-1)=A(k,k-1)P(k-1|k-1)A(k,k-1) T+G(k-1)W(k-1)G(k-1) T
    其中,A(k,k-1)为滤波器k-1时刻到k时刻的滤波器一步转移矩阵,P(k-1|k-1)为k-1时刻的状态估计均方差,P(k|k-1)为k-1时刻到k时刻的一步预测均方差,G(k-1)为滤波器k-1时刻的滤波器噪声系数矩阵,W(k-1)为k-1时刻状态噪声,上标T表示转置;
    2)计算k时刻扩展卡尔曼滤波器滤波增益K(k):
    K(k)=P(k|k-1)H(k) T[H(k)P(k|k-1)H(k) T+R(k)] -1
    其中,K(k)为k时刻的滤波增益,H(k)为k时刻量测矩阵,R(k)为k时刻的量测噪声,上标-1表示求逆;
    3)计算k时刻扩展卡尔曼滤波器状态估计值
    Figure PCTCN2020070025-appb-100032
    Figure PCTCN2020070025-appb-100033
    其中,
    Figure PCTCN2020070025-appb-100034
    为k时刻状态量的估计值,
    Figure PCTCN2020070025-appb-100035
    为k-1到k时刻的状态变量一步预测值,Y(k)为k时刻的量测值,
    Figure PCTCN2020070025-appb-100036
    为k时刻三轴陀螺在X、Y、Z轴上的输出分量;
    4)计算k时刻扩展卡尔曼滤波器估计均方误差P(k|k):
    P(k|k)=[I-K(k)H(k)]P(k|k-1)
    其中,P(k|k)为k时刻估计均方误差,I为单位矩阵;
    5)基于卡尔曼滤波器,通过三轴陀螺的输出Y(k)对状态中的角速度
    Figure PCTCN2020070025-appb-100037
    进行校正。
PCT/CN2020/070025 2019-04-29 2020-01-02 基于角加速度计/陀螺/加速度计的惯性导航解算方法 WO2020220729A1 (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201910354258.4 2019-04-29
CN201910354258.4A CN110207697B (zh) 2019-04-29 2019-04-29 基于角加速度计/陀螺/加速度计的惯性导航解算方法

Publications (1)

Publication Number Publication Date
WO2020220729A1 true WO2020220729A1 (zh) 2020-11-05

Family

ID=67786686

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2020/070025 WO2020220729A1 (zh) 2019-04-29 2020-01-02 基于角加速度计/陀螺/加速度计的惯性导航解算方法

Country Status (2)

Country Link
CN (1) CN110207697B (zh)
WO (1) WO2020220729A1 (zh)

Cited By (18)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112729281A (zh) * 2020-12-18 2021-04-30 无锡卡尔曼导航技术有限公司 一种约束惯性/卫星组合导航静止时航向漂移的方法
CN113295158A (zh) * 2021-05-14 2021-08-24 江苏大学 一种融合惯性数据、地图信息以及行人运动状态的室内定位方法
CN113514064A (zh) * 2021-04-23 2021-10-19 南京航空航天大学 一种鲁棒因子图多源容错导航方法
CN113566849A (zh) * 2021-07-29 2021-10-29 深圳元戎启行科技有限公司 惯性测量单元的安装角度标定方法、装置和计算机设备
CN113624260A (zh) * 2021-08-26 2021-11-09 三一智矿科技有限公司 里程计脉冲当量标定方法及装置、电子设备、存储介质
CN114413934A (zh) * 2022-01-20 2022-04-29 北京经纬恒润科技股份有限公司 一种车辆定位***校正方法和装置
CN114413895A (zh) * 2022-02-24 2022-04-29 中国人民解放军国防科技大学 光纤陀螺旋转惯导联合定位方法、装置、设备及介质
CN114608578A (zh) * 2022-03-07 2022-06-10 北京航空航天大学 一种加权不确定性无人机集群协同导航方法
CN114964228A (zh) * 2022-05-05 2022-08-30 上海机电工程研究所 基于秩亏约束级联滤波的惯性导航协同方法
CN114993296A (zh) * 2022-04-19 2022-09-02 北京自动化控制设备研究所 一种制导炮弹高动态组合导航方法
CN115143993A (zh) * 2022-07-01 2022-10-04 中国人民解放军国防科技大学 基于三轴转台的激光陀螺惯导***g敏感性误差标定方法
CN115265592A (zh) * 2022-07-18 2022-11-01 北京航空航天大学 一种光纤陀螺磁温交联耦合误差的在线补偿方法
CN116046027A (zh) * 2023-03-31 2023-05-02 中国船舶集团有限公司第七〇七研究所 三轴旋转式惯导位置误差无源自主校准方法及***
CN116105731A (zh) * 2023-04-07 2023-05-12 中国人民解放军国防科技大学 稀疏测距条件下的导航方法、装置、计算机设备及介质
CN116136405A (zh) * 2023-04-04 2023-05-19 天津大学 引入磁流体传感器的惯性测量单元的数据处理方法及装置
CN116147577A (zh) * 2023-03-06 2023-05-23 中国人民解放军国防科技大学 基于单轴rins/ldv组合的连续高程测量方法及***
CN116481535A (zh) * 2023-02-02 2023-07-25 中国科学院力学研究所 一种利用惯导数据修正飞行弹道数据的计算方法
CN118067120A (zh) * 2024-04-20 2024-05-24 西安现代控制技术研究所 一种基于主动段激励的高精度惯导航向误差在线估计方法

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110207697B (zh) * 2019-04-29 2023-03-21 南京航空航天大学 基于角加速度计/陀螺/加速度计的惯性导航解算方法
CN113375694A (zh) * 2021-05-25 2021-09-10 南京航空航天大学 一种静基座条件下低成本陀螺全部零偏快速估计方法
CN114279445B (zh) * 2021-12-15 2024-05-24 南京航空航天大学 自旋类飞行器的姿态解算方法

Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2000012966A1 (en) * 1998-07-22 2000-03-09 Litton Systems, Inc. Method for determining gravity in an inertial navigation system
CN102628691A (zh) * 2012-04-09 2012-08-08 北京自动化控制设备研究所 一种完全自主的相对惯性导航方法
CN108562289A (zh) * 2018-06-07 2018-09-21 南京航空航天大学 连续多边几何环境中四旋翼飞行器激光雷达导航方法
CN108592911A (zh) * 2018-03-23 2018-09-28 南京航空航天大学 一种四旋翼飞行器动力学模型/机载传感器组合导航方法
CN108759814A (zh) * 2018-04-13 2018-11-06 南京航空航天大学 一种四旋翼飞行器横滚轴角速度和俯仰轴角速度估计方法
CN108827339A (zh) * 2018-04-10 2018-11-16 南京航空航天大学 一种基于惯性辅助的高效视觉里程计
CN108981708A (zh) * 2018-08-02 2018-12-11 南京航空航天大学 四旋翼扭矩模型/航向陀螺/磁传感器容错组合导航方法
CN108981709A (zh) * 2018-08-02 2018-12-11 南京航空航天大学 基于力矩模型辅助的四旋翼横滚角、俯仰角容错估计方法
CN109612459A (zh) * 2018-11-15 2019-04-12 南京航空航天大学 基于动力学模型的四旋翼飞行器惯性传感器容错导航方法
CN110207697A (zh) * 2019-04-29 2019-09-06 南京航空航天大学 基于角加速度计/陀螺/加速度计的惯性导航解算方法

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2881225B1 (ja) * 1998-03-31 1999-04-12 防衛庁技術研究本部長 慣性航法装置
US8645063B2 (en) * 2010-12-22 2014-02-04 Custom Sensors & Technologies, Inc. Method and system for initial quaternion and attitude estimation
CN102736631B (zh) * 2012-06-11 2015-01-07 北京航空航天大学 一种基于角加速度传感器的多操纵面无人机闭环控制分配方法
CN105203098B (zh) * 2015-10-13 2018-10-02 上海华测导航技术股份有限公司 基于九轴mems传感器的农业机械全姿态角更新方法
CN106767805B (zh) * 2017-01-08 2023-08-18 上海拓攻机器人有限公司 基于mems传感器阵列的高精度惯性量测量方法及测量***
CN106767931A (zh) * 2017-01-24 2017-05-31 北京理工大学 一种基于角加速度给定的无陀螺惯性导航***的验证方法
CN108693372B (zh) * 2018-04-13 2020-07-07 南京航空航天大学 一种四旋翼飞行器的航向轴角速度估计方法

Patent Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2000012966A1 (en) * 1998-07-22 2000-03-09 Litton Systems, Inc. Method for determining gravity in an inertial navigation system
CN102628691A (zh) * 2012-04-09 2012-08-08 北京自动化控制设备研究所 一种完全自主的相对惯性导航方法
CN108592911A (zh) * 2018-03-23 2018-09-28 南京航空航天大学 一种四旋翼飞行器动力学模型/机载传感器组合导航方法
CN108827339A (zh) * 2018-04-10 2018-11-16 南京航空航天大学 一种基于惯性辅助的高效视觉里程计
CN108759814A (zh) * 2018-04-13 2018-11-06 南京航空航天大学 一种四旋翼飞行器横滚轴角速度和俯仰轴角速度估计方法
CN108562289A (zh) * 2018-06-07 2018-09-21 南京航空航天大学 连续多边几何环境中四旋翼飞行器激光雷达导航方法
CN108981708A (zh) * 2018-08-02 2018-12-11 南京航空航天大学 四旋翼扭矩模型/航向陀螺/磁传感器容错组合导航方法
CN108981709A (zh) * 2018-08-02 2018-12-11 南京航空航天大学 基于力矩模型辅助的四旋翼横滚角、俯仰角容错估计方法
CN109612459A (zh) * 2018-11-15 2019-04-12 南京航空航天大学 基于动力学模型的四旋翼飞行器惯性传感器容错导航方法
CN110207697A (zh) * 2019-04-29 2019-09-06 南京航空航天大学 基于角加速度计/陀螺/加速度计的惯性导航解算方法

Cited By (29)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112729281A (zh) * 2020-12-18 2021-04-30 无锡卡尔曼导航技术有限公司 一种约束惯性/卫星组合导航静止时航向漂移的方法
CN113514064A (zh) * 2021-04-23 2021-10-19 南京航空航天大学 一种鲁棒因子图多源容错导航方法
CN113514064B (zh) * 2021-04-23 2024-01-30 南京航空航天大学 一种鲁棒因子图多源容错导航方法
CN113295158A (zh) * 2021-05-14 2021-08-24 江苏大学 一种融合惯性数据、地图信息以及行人运动状态的室内定位方法
CN113295158B (zh) * 2021-05-14 2024-05-14 江苏大学 一种融合惯性数据、地图信息以及行人运动状态的室内定位方法
CN113566849B (zh) * 2021-07-29 2024-03-05 深圳元戎启行科技有限公司 惯性测量单元的安装角度标定方法、装置和计算机设备
CN113566849A (zh) * 2021-07-29 2021-10-29 深圳元戎启行科技有限公司 惯性测量单元的安装角度标定方法、装置和计算机设备
CN113624260A (zh) * 2021-08-26 2021-11-09 三一智矿科技有限公司 里程计脉冲当量标定方法及装置、电子设备、存储介质
CN113624260B (zh) * 2021-08-26 2024-02-27 三一智矿科技有限公司 里程计脉冲当量标定方法及装置、电子设备、存储介质
CN114413934A (zh) * 2022-01-20 2022-04-29 北京经纬恒润科技股份有限公司 一种车辆定位***校正方法和装置
CN114413934B (zh) * 2022-01-20 2024-01-26 北京经纬恒润科技股份有限公司 一种车辆定位***校正方法和装置
CN114413895A (zh) * 2022-02-24 2022-04-29 中国人民解放军国防科技大学 光纤陀螺旋转惯导联合定位方法、装置、设备及介质
CN114413895B (zh) * 2022-02-24 2023-10-13 中国人民解放军国防科技大学 光纤陀螺旋转惯导联合定位方法、装置、设备及介质
CN114608578A (zh) * 2022-03-07 2022-06-10 北京航空航天大学 一种加权不确定性无人机集群协同导航方法
CN114608578B (zh) * 2022-03-07 2024-04-09 北京航空航天大学 一种加权不确定性无人机集群协同导航方法
CN114993296A (zh) * 2022-04-19 2022-09-02 北京自动化控制设备研究所 一种制导炮弹高动态组合导航方法
CN114993296B (zh) * 2022-04-19 2024-03-15 北京自动化控制设备研究所 一种制导炮弹高动态组合导航方法
CN114964228A (zh) * 2022-05-05 2022-08-30 上海机电工程研究所 基于秩亏约束级联滤波的惯性导航协同方法
CN115143993B (zh) * 2022-07-01 2024-03-29 中国人民解放军国防科技大学 基于三轴转台的激光陀螺惯导***g敏感性误差标定方法
CN115143993A (zh) * 2022-07-01 2022-10-04 中国人民解放军国防科技大学 基于三轴转台的激光陀螺惯导***g敏感性误差标定方法
CN115265592A (zh) * 2022-07-18 2022-11-01 北京航空航天大学 一种光纤陀螺磁温交联耦合误差的在线补偿方法
CN115265592B (zh) * 2022-07-18 2024-04-09 北京航空航天大学 一种光纤陀螺磁温交联耦合误差的在线补偿方法
CN116481535A (zh) * 2023-02-02 2023-07-25 中国科学院力学研究所 一种利用惯导数据修正飞行弹道数据的计算方法
CN116147577B (zh) * 2023-03-06 2024-05-03 中国人民解放军国防科技大学 基于单轴rins/ldv组合的连续高程测量方法及***
CN116147577A (zh) * 2023-03-06 2023-05-23 中国人民解放军国防科技大学 基于单轴rins/ldv组合的连续高程测量方法及***
CN116046027A (zh) * 2023-03-31 2023-05-02 中国船舶集团有限公司第七〇七研究所 三轴旋转式惯导位置误差无源自主校准方法及***
CN116136405A (zh) * 2023-04-04 2023-05-19 天津大学 引入磁流体传感器的惯性测量单元的数据处理方法及装置
CN116105731A (zh) * 2023-04-07 2023-05-12 中国人民解放军国防科技大学 稀疏测距条件下的导航方法、装置、计算机设备及介质
CN118067120A (zh) * 2024-04-20 2024-05-24 西安现代控制技术研究所 一种基于主动段激励的高精度惯导航向误差在线估计方法

Also Published As

Publication number Publication date
CN110207697B (zh) 2023-03-21
CN110207697A (zh) 2019-09-06

Similar Documents

Publication Publication Date Title
WO2020220729A1 (zh) 基于角加速度计/陀螺/加速度计的惯性导航解算方法
CN109813311B (zh) 一种无人机编队协同导航方法
CN109596018B (zh) 基于磁测滚转角速率信息的旋转弹飞行姿态高精度估计方法
CN106643737B (zh) 风力干扰环境下四旋翼飞行器姿态解算方法
CN104698485B (zh) 基于bd、gps及mems的组合导航***及导航方法
CN102809377A (zh) 飞行器惯性/气动模型组合导航方法
CN111207745B (zh) 一种适用于大机动无人机垂直陀螺仪的惯性测量方法
CN108007477B (zh) 一种基于正反向滤波的惯性行人定位***误差抑制方法
CN110954102B (zh) 用于机器人定位的磁力计辅助惯性导航***及方法
CN111024070A (zh) 一种基于航向自观测的惯性足绑式行人定位方法
CN105865455A (zh) 一种利用gps与加速度计计算飞行器姿态角的方法
CN113340298B (zh) 一种惯导和双天线gnss外参标定方法
CN108534783A (zh) 一种基于北斗导航技术的飞行器导航方法
CN108592911A (zh) 一种四旋翼飞行器动力学模型/机载传感器组合导航方法
CN108592943A (zh) 一种基于opreq方法的惯性系粗对准计算方法
CN112556724A (zh) 动态环境下的微型飞行器低成本导航***初始粗对准方法
CN113959462B (zh) 一种基于四元数的惯性导航***自对准方法
CN111207773A (zh) 一种用于仿生偏振光导航的姿态无约束优化求解方法
CN107063248A (zh) 基于旋翼转速的动力学模型辅助惯导的导航方法
CN111141285B (zh) 一种航空重力测量装置
CN113932803B (zh) 适用于高动态飞行器的惯性/地磁/卫星组合导航***
Du et al. A low-cost attitude estimation system for UAV application
Yang et al. Model-free integrated navigation of small fixed-wing UAVs full state estimation in wind disturbance
WO2021012635A1 (zh) 一种基于陀螺仪信息的惯性导航方法
Cao et al. Research of attitude estimation of UAV based on information fusion of complementary filter

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 20797981

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20797981

Country of ref document: EP

Kind code of ref document: A1