CN101706287A - Rotating strapdown system on-site proving method based on digital high-passing filtering - Google Patents

Rotating strapdown system on-site proving method based on digital high-passing filtering Download PDF

Info

Publication number
CN101706287A
CN101706287A CN200910073242A CN200910073242A CN101706287A CN 101706287 A CN101706287 A CN 101706287A CN 200910073242 A CN200910073242 A CN 200910073242A CN 200910073242 A CN200910073242 A CN 200910073242A CN 101706287 A CN101706287 A CN 101706287A
Authority
CN
China
Prior art keywords
imu
omega
state
matrix
order
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
CN200910073242A
Other languages
Chinese (zh)
Other versions
CN101706287B (en
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.)
Harbin Engineering University
Original Assignee
Harbin Engineering 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 Harbin Engineering University filed Critical Harbin Engineering University
Priority to CN2009100732422A priority Critical patent/CN101706287B/en
Publication of CN101706287A publication Critical patent/CN101706287A/en
Application granted granted Critical
Publication of CN101706287B publication Critical patent/CN101706287B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Navigation (AREA)

Abstract

The invention provides a rotating strapdown system on-site proving method based on digital high-passing filtering, which comprises the steps: (1) determining an initial position parameter of a carrier by a GPS; (2) collecting data output by an optical fiber gyroscope and an acceleration meter, and processing the data; (3) rotating and stopping at four positions of a single axis of an inertia measurement unit; (4) analyzing biased observability degree of an inertia device by using a spectral condition number method; (5) filtering the Schuler period contained in velocity information of a navigation system by adopting an IIR high-passing digital filter; and (6) taking filtered velocity information as observed quantity, and estimating the deviation of the inertia device by adopting a Kalman filtering technique. When the carrier is in anchoring state, the adoption of the method can obtain higher on-site proving precision.

Description

A kind of rotating strapdown system field calibration method based on Digital High Pass Filter
Technical field
What the present invention relates to is a kind of measuring method, in particular a kind of rotating strapdown system field calibration method based on Digital High Pass Filter.
Background technology
Strapdown inertial navigation system (SINS) is connected in inertance elements such as gyroscope, accelerometer on the carrier, according to Newton mechanics law, by the information of these inertance element collections is handled, obtain the complete independent navigation equipment of the full navigation information such as attitude, speed, position, acceleration, angular velocity and angular acceleration of carrier.Because SINS relies on its own inertial element fully, do not rely on any external information to measure navigational parameter, therefore, it has good concealment, be not subjected to the weather condition restriction, advantage such as interference-free is a kind of complete autonomous type, round-the-clock navigational system, has been widely used in fields such as Aeronautics and Astronautics, navigation.According to the ultimate principle of SINS, the existence that SINS inertia device in navigation procedure often is worth deviation is the principal element that causes the inertial navigation system navigation accuracy to be difficult to improve.How effectively limiting the inertial navigation error, to disperse, improve the inertial navigation system precision be the very important problem in inertial navigation field.
The error of inertial navigation system suppresses, is not to depend on outside assisting error state is estimated, but the propagation law of research inertial navigation error under the special exercise condition, and limit error according to this rule and disperse, improve the method for navigation accuracy.Rotating inhibition is most typical error inhibition method: by around an axle or a plurality of rotator inertia measuring units (IMU), navigation error is modulated, reach the purpose that navigation accuracy is dispersed, improved to the control navigation error.
Calibration technique is exactly a kind of method that improves the actual service precision of inertial navigation from the software aspect.Calibration technique also is a kind of Error Compensation Technology in essence.So-called Error Compensation Technology is exactly to set up the error mathematic model of inertance element and inertial navigation system, determines model coefficient by certain test, and then eliminates error by software algorithm.Inertance element and inertial navigation system must be by demarcating to determine basic error mathematic model parameter, to guarantee the operate as normal of element and system before dispatching from the factory.And the error compensation under the abominable dynamic environment of research, inertial navigation system of inertance element high-order error term all is to carry out on the basis of demarcating, and we can say that staking-out work is the basis of whole Error Compensation Technology.The Inertial Measurement Unit error (drift of inertance element and scale factor error) of starting shooting one by one is different, and along with the time increases the IMU output error and drifts about in time, realizes that on-the-spot on-line proving is significant for improving system accuracy.
Open report related to the present invention in the CNKI storehouse has: 1, " horizontal initial alignment error is to the influence of rotation IMU accuracy of navigation systems ", this article has mainly been told about the influence of horizontal misalignment to the rotation strapdown inertial navitation system (SINS), tell about the unidirectional rotation of IMU single shaft in its Chinese, but do not mentioned the content of on-site proving.2, " application of rotation IMU in fiber strapdown boat appearance system ", this article has mainly been introduced single shaft, twin shaft rotation mode, and proves in theory.3, " the six positions rotation on-site proving new method of optical fibre gyro IMU ", this article has mainly been introduced at the scene of using optical fibre gyro IMU has been carried out the rotation of ten secondaries on six position, set up 42 non-linear input and output equations according to the error model of optical fibre gyro IMU then, disappear mutually with the symmetric position error by the rotation integration, eliminate the nonlinear terms in the equation, finally solve gyro constant multiplier, gyroscope constant value drift, gyro misalignment and accelerometer and often be worth 15 error coefficients such as biasing.
Summary of the invention
The object of the present invention is to provide a kind of carrier of working as to be under the moored condition, can obtain a kind of rotating strapdown system field calibration method of higher on-site proving precision based on Digital High Pass Filter.
The object of the present invention is achieved like this: may further comprise the steps:
(1) utilizes global position system GPS to determine the initial position parameters of carrier, they are bound to navigational computer;
(2) fiber optic gyro strapdown inertial navigation system carries out gathering after the preheating data of fibre optic gyroscope and quartz accelerometer output;
(3) IMU adopts 8 commentaries on classics to stop the transposition scheme that order is a swing circle;
(4) utilize the spectral condition method of counting to ask for the observability degree of inertia device deviation in the IMU four-position rotation and stop process;
(5) designing the infinite impulse response digital high-pass filter, is that the carrier horizontal velocity that calculates is down carried out the high-pass filtering processing with navigating, and the Schuler period under filtering is navigated and is in the bearer rate keeps carrier owing to the velocity deviation of waving and swinging the generation of moving;
Model is floated in estimating when (6) setting up the carrier moored condition according to the moving pedestal error equation of inertial navigation system, directly as observed quantity, utilizes Kalman Filter Technology to realize the on-site proving of rotation strapdown inertial navitation system (SINS) with the speed that obtains after the high-pass filtering.
The present invention can also comprise:
1, to adopt 8 commentaries on classics to stop order be that the transposition scheme of a swing circle is for described IMU:
Order 1, IMU clockwise rotates 180 ° of in-position C, stand-by time T from the A point tOrder 2, IMU clockwise rotates 90 ° of in-position D, stand-by time T from the C point tOrder 3, IMU rotates counterclockwise 180 ° of in-position B, stand-by time T from the D point tOrder 4, IMU rotates counterclockwise 90 ° of in-position A, stand-by time T from the B point tOrder 5, IMU rotates counterclockwise 180 ° of in-position C, stand-by time T from the A point tOrder 6, IMU rotates counterclockwise 90 ° of in-position B, stand-by time T from the C point tOrder 7, IMU clockwise rotates 180 ° of in-position D, stand-by time T from the B point tOrder 8, IMU clockwise rotates 90 ° of in-position A, stand-by time T from the D point tIMU rotates sequential loop according to this to carry out; IMU pause point p3, p8 and p4, p7 on the horizontal east orientation axle are symmetrical in the rotating shaft center; Pause point p1 on the north orientation axle, p5 and p2, p6 are symmetrical in the rotating shaft center; It is that carry out at 180 ° or 90 ° of intervals that four-position rotation and stop scheme remains rotational angle.
2, the described method of utilizing the spectral condition method of counting to ask for the observability degree of inertia device deviation in the IMU four-position rotation and stop process is:
Find the solution system of linear equations
AX=b,b∈C n
If A ∈ is C N * n, || || be a kind of operator norm,
cond ( A ) = | | A | | | | A - 1 | | , det A ≠ 0 ∞ , det A = 0
Claim cond (A) be matrix A about operator norm || || conditional number, commonly used is about the p-norm || || pConditional number, note is made cond p(A), cond 2(A) be the spectral condition number,
At discrete time-varying system:
X k + 1 = F k X k Z k = H X k
Bring system state equation into observation equation and obtain a set of equations:
Z 0 = H X 0 Z 1 = H F 0 X 0 . . . Z k = H Π i = 0 k F i X 0
Note
O k = H HF 0 . . . H Π i = 0 k F i T
Then
O kX 0=Z
For stational system F kBe constant, O kBe exactly the ornamental matrix, time-varying system is observed on sampled point and is obtained discrete time-varying system, F kBe exactly the state-transition matrix Φ (t in the sampling period k+ T, t k),
O k=[HHΦ(t 1,t 0)…HΦ(t k,t 0)] T
State is the n dimension, an observed quantity Z kBe r dimension (r<n), the order of observation battle array H is r, carries out k time at least and observes (kr 〉=n), obtain X 0, find the solution state X according to least square method 0,
X 0 = ( O k T O k ) - 1 O k T Z
Figure G2009100732422D0000044
Be the observation battle array, because Be normal matrix, count cond by calculating spectral condition 2(M), analyze stability of solution, and
cond 2 ( M ) = max λ ∈ λ ( M ) | λ | min λ ∈ λ ( M ) | λ |
In the formula, λ is the eigenwert of matrix M, eigenwert and the proper vector of further analysis matrix M, so that determine that the observability degree of which state is better actually, the observability degree of which state is poor, but with M diagonalization at the tenth of the twelve Earthly Branches, note U TMU=Λ, wherein Λ=diag (λ 1, λ 2... λ n), the observability degree S of state X then:
S=abs(U[λ 1,λ 2,...,λ n] T)
Calculate eigenwert and the proper vector of the observability matrix M of system, determine the observability degree of each state.
3, the described method of utilizing Kalman Filter Technology to realize the on-site proving of rotation strapdown inertial navitation system (SINS) is:
1) set up the state equation of Kalman filtering:
The state error of rotation strapdown inertial navitation system (SINS) is described with linear first-order differential equation:
X · ( t ) = A ( t ) X ( t ) + B ( t ) W ( t )
The state vector of etching system when wherein X (t) is t; A (t) and B (t) are respectively the state matrix and the noise matrix of system; W (t) is the system noise vector;
The state vector of system is:
Figure G2009100732422D0000052
The white noise vector of system is:
W=[a x?a yxyz?0?0?0?0?0?0?0?0] T
δ V wherein e, δ V nThe velocity error of representing east orientation, north orientation respectively; Be respectively IMU coordinate system ox s, oy sAxis accelerometer zero partially; ε x, ε y, ε zBe respectively IMU coordinate system ox s, oy s, oz sThe constant value drift of axle gyro; a x, a yBe respectively IMU coordinate system ox s, oy sThe white noise error of axis accelerometer; δ K Gx, δ K Gy, δ K GzBe respectively IMU coordinate system ox s, oy s, oz sThe constant multiplier error of axle gyro; ω x, ω y, ω zBe respectively IMU coordinate system ox s, oy s, oz sThe white noise error of axle gyro;
The state-transition matrix of system is:
A = F 2 × 2 1 f 2 × 3 T ~ 2 × 2 O 2 × 6 F 3 × 2 2 F 3 × 3 3 O 3 × 2 T 3 × 6 O 8 × 2 O 8 × 3 O 8 × 2 O 8 × 6
F 2 × 2 1 = V N R n tan L 2 ω ie sin L + V E R n tan L - ( 2 ω ie sin L + 2 V E R n tan L ) 0
F 3 × 2 2 = 0 - 1 R m 1 R n 0 tan L R n 0
F 3 × 3 3 = 0 ω ie sin L + V E tan L R n - ( ω ie cos L + V E R n ) - ( ω ie sin L + V E tan L R n ) 0 - V N R m ω ie cos L + V E R n V N R m 0
f 2 × 3 = 0 - f U f N f U 0 f E
T ~ 2 × 2 = T 11 T 12 T 21 T 22
T 3 × 6 = - T 11 - T 12 - T 13 - T 11 ω x - T 12 ω y - T 13 ω z - T 21 - T 22 - T 23 - T 21 ω x - T 22 ω y - T 23 ω z - T 31 - T 32 - T 33 - T 31 ω x - T 32 ω y - T 33 ω z
V E, V NThe speed of representing east orientation, north orientation respectively; ω x, ω y, ω zThree input angular velocities representing gyro respectively; ω IeThe expression rotational-angular velocity of the earth; R m, R nRepresent earth meridian, fourth of the twelve Earthly Branches radius-of-curvature at the tenth of the twelve Earthly Branches respectively; L represents local latitude; f E, f N, f UBe expressed as respectively navigation coordinate system down east orientation, north orientation, day to specific force;
2) set up the measurement equation of Kalman filtering:
The measurement equation of describing the rotation strapdown inertial navitation system (SINS) with linear first-order differential equation is as follows:
Z(t)=H(t)X(t)+V(t)
Wherein: the measurement vector of etching system during Z (t) expression t; The measurement matrix of H (t) expression system; The measurement noise of V (t) expression system;
The system measurements matrix is:
H = 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
Amount is measured as the speed that obtains after the high-pass filtering:
Z = V ~ E V ~ N T .
Four positions that the present invention is fixing with the relative carrier azimuth axis of Inertial Measurement Unit are changeed and are stopped, utilize the moving observability degree that improves systematic parameter of commentaries on classics stoppage in transit of IMU, the design Kalman filter is introduced each calibrated error item that the high precision oracle encourages IMU respectively, estimate and compensation IMU output error, finish the on-site proving of system.
The present invention's advantage compared with prior art is: the present invention has broken the on-site proving that traditional scaling method is not suitable for system, proposing a kind of IMU of utilization commentaries on classics stops improving the systematic parameter observability degree and adopts Kalman Filter Technology the inertia device deviation to be carried out the scheme of on-line proving, this method inertia device often can be worth deviation and gyroscope constant multiplier error is carried out scene estimation and compensation, can improve system alignment and navigation accuracy effectively.
The effect useful to the present invention is described as follows:
Under Visual C++ simulated conditions, this method is carried out emulation experiment:
Setting the gyroscope constant value drift is 0.01 °/h, and the accelerometer zero drift is 10 -4G; System's initial alignment error is 0.1 °, 0.1 °, 0.5 °; Carrier is done the three-axis swinging motion with sinusoidal rule around pitch axis, axis of roll and course axle, and its mathematical model is:
θ = θ m sin ( ω θ t + φ θ ) γ = γ m sin ( ω γ t + φ γ ) ψ = ψ m sin ( ω ψ t + φ ψ ) + k
Wherein: θ, γ, ψ represent the angle variables of waving of pitch angle, roll angle and course angle respectively; θ m, γ m, ψ mThe angle amplitude is waved in expression accordingly respectively; ω θ, ω γ, ω ψRepresent corresponding angle of oscillation frequency respectively; φ θ, φ γ, φ ψRepresent corresponding initial phase respectively; ω i=2 τ/T i, i=θ, γ, ψ, T iRepresent corresponding rolling period, k is the angle, initial heading.Get during emulation: θ m=6 °, γ m=12 °, ψ m=10 °, T θ=8s, T γ=10s, T ψ=6s, k=0 °.
The swaying of carrier, surging and hang down and swing the acceleration that causes and be:
Figure G2009100732422D0000072
In the formula, i=x, y, z be geographic coordinate system east orientation, north orientation, day to.
Figure G2009100732422D0000073
Figure G2009100732422D0000075
Figure G2009100732422D0000077
Figure G2009100732422D0000078
Figure G2009100732422D0000079
For going up, [0,2 π] obey equally distributed random phase.
Carrier initial position: 45.7796 ° of north latitude, 126.6705 ° of east longitudes;
The initial attitude error angle: three initial attitude error angles are zero;
Equatorial radius: R e=6378393.0m;
Ellipsoid degree: e=3.367e-3;
The earth surface acceleration of gravity that can get by universal gravitation: g 0=9.78049;
Rotational-angular velocity of the earth (radian per second): 7.2921158e-5;
The gyroscope constant value drift: 0.01 degree/hour;
The gyroscope random walk: 0.001 degree/
Figure G2009100732422D0000081
Gyroscope constant multiplier error: 1000ppm;
Accelerometer bias: 10 -4g 0
Accelerometer noise: 10 -6g 0
Constant: π=3.1415926;
The mathematical model parameter of IMU four-position rotation and stop scheme:
The dead time of four positions: T t=5min;
The time that consumes when rotating 180 ° and 90 °: T z=12s;
Rotate in the process of 180 ° and 90 °, the acceleration and deceleration time in each transposition respectively is 4s.
Description of drawings
Fig. 1 is a kind of rotating strapdown system field calibration method process flow diagram based on Digital High Pass Filter of the present invention;
Fig. 2 is in the IMU four-position rotation and stop process of the present invention, the relative position relation of IMU coordinate system and carrier coordinate system;
Fig. 3 is the gyroscope constant value drift of estimation of the present invention;
Fig. 4 is the accelerometer bias on the horizontal direction of estimation of the present invention;
Fig. 5 is the gyroscope constant multiplier error of estimation of the present invention.
Embodiment
For example the present invention is done description in more detail below in conjunction with accompanying drawing:
(1) utilizes global position system GPS to determine the initial position parameters of carrier, they are bound to navigational computer;
(2) fiber optic gyro strapdown inertial navigation system carries out gathering after the preheating data of fibre optic gyroscope and quartz accelerometer output;
(3) IMU adopts 8 commentaries on classics to stop the transposition scheme that order is a swing circle;
Order 1, IMU clockwise rotates 180 ° of in-position C, stand-by time T from the A point tOrder 2, IMU clockwise rotates 90 ° of in-position D, stand-by time T from the C point tOrder 3, IMU rotates counterclockwise 180 ° of in-position B, stand-by time T from the D point tOrder 4, IMU rotates counterclockwise 90 ° of in-position A, stand-by time T from the B point tOrder 5, IMU rotates counterclockwise 180 ° of in-position C, stand-by time T from the A point tOrder 6, IMU rotates counterclockwise 90 ° of in-position B, stand-by time T from the C point tOrder 7, IMU clockwise rotates 180 ° of in-position D, stand-by time T from the B point tOrder 8, IMU clockwise rotates 90 ° of in-position A, stand-by time T from the D point tIMU rotates sequential loop according to this to carry out.Positive and negative average in order effectively the inertia device deviation on the horizontal direction to be carried out on symmetric position, IMU pause point p3, p8 and p4, p7 on the horizontal east orientation axle of definition are symmetrical in the rotating shaft center; Pause point p1 on the north orientation axle, p5 and p2, p6 are symmetrical in the rotating shaft center.It is that carry out at 180 ° or 90 ° of intervals that improved four-position rotation and stop scheme remains rotational angle.
(4) utilize the spectral condition method of counting to ask for the observability degree of inertia device deviation in the IMU four-position rotation and stop process;
Consider to find the solution system of linear equations
AX=b,b∈C n (1)
If A ∈ is C N * n, || || be a kind of operator norm,
cond ( A ) = | | A | | | | A - 1 | | , det A ≠ 0 ∞ , det A = 0 - - - ( 2 )
Claim cond (A) be matrix A (about inverting or finding the solution system of linear equations) about operator norm || || conditional number.Commonly used is about the p-norm || || pConditional number, can remember and make cond p(A).Especially, claim cond 2(A) be the spectral condition number.
At discrete time-varying system:
X k + 1 = F k X k Z k = H X k - - - ( 3 )
Can be by a group observations Z=[Z 0, Z 1... .Z k] TObtain original state X 0Be the considerable essence of system.Bring system state equation into observation equation and obtain a set of equations:
Z 0 = H X 0 Z 1 = H F 0 X 0 . . . Z k = H Π i = 0 k F i X 0 - - - ( 4 )
Note
O k = H HF 0 . . . H Π i = 0 k F i T - - - ( 5 )
Then have
O kX 0=Z (6)
For stational system F kBe constant, O kIt is exactly the ornamental matrix.Time-varying system is observed on sampled point and is obtained discrete time-varying system, F kBe exactly the state-transition matrix Φ (t in the sampling period k+ T, t k), therefore have
O k=[H?HΦ(t 1,t 0)…HΦ(t k,t 0)] T (7)
If state is the n dimension, an observed quantity Z kBe r dimension (r<n), the order of observation battle array H is r, carries out k time at least and observes (kr 〉=n), just can obtain X 0Find the solution state X according to least square method 0
X 0 = ( O k T O k ) - 1 O k T Z - - - ( 8 )
Note
Figure G2009100732422D0000105
Be the observation battle array.Because
Figure G2009100732422D0000106
Be normal matrix, then can count cond by calculating spectral condition 2(M), analyze stability of solution, and
cond 2 ( M ) = max λ ∈ λ ( M ) | λ | min λ ∈ λ ( M ) | λ | - - - ( 9 )
In the formula, λ is the eigenwert of matrix M.Eigenwert and the proper vector of further analysis matrix M, so that determine that the observability degree of which state is better actually, the observability degree of which state is poor.But with M diagonalization at the tenth of the twelve Earthly Branches, note U TMU=Λ, wherein Λ=diag (λ 1, λ 2... λ n), the observability degree S of state X then:
S=abs(U[λ 1,λ 2,...,λ n] T) (10)
Calculate eigenwert and the proper vector of the observability matrix M of system, just can determine the observability degree of each state.
(5) design infinite impulse response digital high-pass filter (IIR), the carrier horizontal velocity that navigation system calculates is down carried out the high-pass filtering processing, Schuler period under filtering is navigated and is in the bearer rate keeps carrier owing to wave and swing the velocity deviation of the generation of moving;
Model is floated in estimating when (6) setting up the carrier moored condition according to the moving pedestal error equation of inertial navigation system, with the speed that obtains after the high-pass filtering directly as observed quantity.Utilize Kalman Filter Technology to realize the on-site proving of rotation strapdown inertial navitation system (SINS);
Foundation is the Kalman filter model of observed quantity with the horizontal velocity under being through the navigation after the high-pass filtering;
1) set up the state equation of Kalman filtering:
The state error of rotation strapdown inertial navitation system (SINS) is described with linear first-order differential equation:
X · ( t ) = A ( t ) X ( t ) + B ( t ) W ( t ) - - - ( 11 )
The state vector of etching system when wherein X (t) is t; A (t) and B (t) are respectively the state matrix and the noise matrix of system; W (t) is the system noise vector;
The state vector of system is:
Figure G2009100732422D0000112
The white noise vector of system is:
W=[a x?a yxyz?0?0?0?0?0?0?0?0] T (13)
δ V wherein e, δ V nThe velocity error of representing east orientation, north orientation respectively;
Figure G2009100732422D0000113
Be respectively IMU coordinate system ox s, oy sAxis accelerometer zero partially; ε x, ε y, ε zBe respectively IMU coordinate system ox s, oy s, oz sThe constant value drift of axle gyro; a x, a yBe respectively IMU coordinate system ox s, oy sThe white noise error of axis accelerometer; δ K Gx, δ K Gy, δ K GzBe respectively IMU coordinate system ox s, oy s, oz sThe constant multiplier error of axle gyro; ω x, ω y, ω zBe respectively IMU coordinate system ox s, oy s, oz sThe white noise error of axle gyro;
The state-transition matrix of system is:
A = F 2 × 2 1 f 2 × 3 T ~ 2 × 2 O 2 × 6 F 3 × 2 2 F 3 × 3 3 O 3 × 2 T 3 × 6 O 8 × 2 O 8 × 3 O 8 × 2 O 8 × 6 - - - ( 14 )
F 2 × 2 1 = V N R n tan L 2 ω ie sin L + V E R n tan L - ( 2 ω ie sin L + 2 V E R n tan L ) 0 - - - ( 15 )
F 3 × 2 2 = 0 - 1 R m 1 R n 0 tan L R n 0 - - - ( 16 )
F 3 × 3 3 = 0 ω ie sin L + V E tan L R n - ( ω ie cos L + V E R n ) - ( ω ie sin L + V E tan L R n ) 0 - V N R m ω ie cos L + V E R n V N R m 0 - - - ( 17 )
f 2 × 3 = 0 - f U f N f U 0 f E - - - ( 18 )
T ~ 2 × 2 = T 11 T 12 T 21 T 22 - - - ( 19 )
T 3 × 6 = - T 11 - T 12 - T 13 - T 11 ω x - T 12 ω y - T 13 ω z - T 21 - T 22 - T 23 - T 21 ω x - T 22 ω y - T 23 ω z - T 31 - T 32 - T 33 - T 31 ω x - T 32 ω y - T 33 ω z - - - ( 20 )
V E, V NThe speed of representing east orientation, north orientation respectively; ω x, ω y, ω zThree input angular velocities representing gyro respectively; ω IeThe expression rotational-angular velocity of the earth; R m, R nRepresent earth meridian, fourth of the twelve Earthly Branches radius-of-curvature at the tenth of the twelve Earthly Branches respectively; L represents local latitude; f E, f N, f UBe expressed as respectively navigation coordinate system down east orientation, north orientation, day to specific force.
2) set up the measurement equation of Kalman filtering:
The measurement equation of describing the rotation strapdown inertial navitation system (SINS) with linear first-order differential equation is as follows:
Z(t)=H(t)X(t)+V(t) (21)
Wherein: the measurement vector of etching system during Z (t) expression t; The measurement matrix of H (t) expression system; The measurement noise of V (t) expression system;
The system measurements matrix is:
H = 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 - - - ( 22 )
Amount is measured as the speed that obtains after the high-pass filtering:
Z = V ~ E V ~ N T - - - ( 23 )

Claims (4)

1. rotating strapdown system field calibration method based on Digital High Pass Filter is characterized in that may further comprise the steps:
(1) utilizes global position system GPS to determine the initial position parameters of carrier, they are bound to navigational computer;
(2) fiber optic gyro strapdown inertial navigation system carries out gathering after the preheating data of fibre optic gyroscope and quartz accelerometer output;
(3) IMU adopts 8 commentaries on classics to stop the transposition scheme that order is a swing circle;
(4) utilize the spectral condition method of counting to ask for the observability degree of inertia device deviation in the IMU four-position rotation and stop process;
(5) designing the infinite impulse response digital high-pass filter, is that the carrier horizontal velocity that calculates is down carried out the high-pass filtering processing with navigating, and the Schuler period under filtering is navigated and is in the bearer rate keeps carrier owing to the velocity deviation of waving and swinging the generation of moving;
Model is floated in estimating when (6) setting up the carrier moored condition according to the moving pedestal error equation of inertial navigation system, directly as observed quantity, utilizes Kalman Filter Technology to realize the on-site proving of rotation strapdown inertial navitation system (SINS) with the speed that obtains after the high-pass filtering.
2. a kind of rotating strapdown system field calibration method based on Digital High Pass Filter according to claim 1 is characterized in that it is that the transposition scheme of a swing circle is that described IMU adopts 8 commentaries on classics to stop order:
Order 1, IMU clockwise rotates 180 ° of in-position C, stand-by time T from the A point tOrder 2, IMU clockwise rotates 90 ° of in-position D, stand-by time T from the C point tOrder 3, IMU rotates counterclockwise 180 ° of in-position B, stand-by time T from the D point tOrder 4, IMU rotates counterclockwise 90 ° of in-position A, stand-by time T from the B point tOrder 5, IMU rotates counterclockwise 180 ° of in-position C, stand-by time T from the A point tOrder 6, IMU rotates counterclockwise 90 ° of in-position B, stand-by time T from the C point tOrder 7, IMU clockwise rotates 180 ° of in-position D, stand-by time T from the B point tOrder 8, IMU clockwise rotates 90 ° of in-position A, stand-by time T from the D point tIMU rotates sequential loop according to this to carry out; IMU pause point p3, p8 and p4, p7 on the horizontal east orientation axle are symmetrical in the rotating shaft center; Pause point p1 on the north orientation axle, p5 and p2, p6 are symmetrical in the rotating shaft center; It is that carry out at 180 ° or 90 ° of intervals that four-position rotation and stop scheme remains rotational angle.
3. a kind of rotating strapdown system field calibration method based on Digital High Pass Filter according to claim 2 is characterized in that the described method of utilizing the spectral condition method of counting to ask for the observability degree of inertia device deviation in the IMU four-position rotation and stop process is:
Find the solution system of linear equations
AX=b,b∈C n
If A ∈ is C N * n, || || be a kind of operator norm,
cond ( A ) = | | A | | | | A - 1 | | , det A ≠ 0 ∞ , det A = 0
Claim cond (A) be matrix A about operator norm || || conditional number, commonly used is about the p-norm || || pConditional number, note is made cond p(A), cond 2(A) be the spectral condition number,
At discrete time-varying system:
X k + 1 = F k X k Z k = HX k
Bring system state equation into observation equation and obtain a set of equations:
Z 0 = HX 0 Z 1 = HF 0 X 0 . . . Z k = H Π i = 0 k F i X 0
Note
O k = H HF 0 . . . H Π i = 0 k F i T
Then
O kX 0=Z
For stational system F kBe constant, O kBe exactly the ornamental matrix, time-varying system is observed on sampled point and is obtained discrete time-varying system, F kBe exactly the state-transition matrix Φ (t in the sampling period k+ T, t k),
O k=[H?HΦ(t 1,t 0)…HΦ(t k,t 0)] T
State is the n dimension, an observed quantity Z kBe r dimension (r<n), the order of observation battle array H is r, carries out k time at least and observes (kr 〉=n), obtain X 0, find the solution state X according to least square method 0,
X 0 = ( O k T O k ) - 1 O k T Z
Figure F2009100732422C0000032
Be the observation battle array, because
Figure F2009100732422C0000033
Be normal matrix, count cond by calculating spectral condition 2(M), analyze stability of solution, and
cond 2 ( M ) = max λ ∈ λ ( M ) | λ | min λ ∈ λ ( M ) | λ |
In the formula, λ is the eigenwert of matrix M, eigenwert and the proper vector of further analysis matrix M, so that determine that the observability degree of which state is better actually, the observability degree of which state is poor, but with M diagonalization at the tenth of the twelve Earthly Branches, note U TMU=Λ, wherein Λ=diag (λ 1, λ 2... λ n), the observability degree S of state X then:
S=abs(U[λ 1,λ 2,...,λ n] T)
Calculate eigenwert and the proper vector of the observability matrix M of system, determine the observability degree of each state.
4. a kind of rotating strapdown system field calibration method based on Digital High Pass Filter according to claim 3 is characterized in that the described method of utilizing Kalman Filter Technology to realize the on-site proving of rotation strapdown inertial navitation system (SINS) is:
1) set up the state equation of Kalman filtering:
The state error of rotation strapdown inertial navitation system (SINS) is described with linear first-order differential equation:
X · ( t ) = A ( t ) X ( t ) + B ( t ) W ( t )
The state vector of etching system when wherein X (t) is t; A (t) and B (t) are respectively the state matrix and the noise matrix of system; W (t) is the system noise vector;
The state vector of system is:
The white noise vector of system is:
W=[a x?a yxyz?0?0?0?0?0?0?0?0] T
δ V wherein e, δ V nThe velocity error of representing east orientation, north orientation respectively;
Figure F2009100732422C0000041
Be respectively IMU coordinate system ox s, oy sAxis accelerometer zero partially; ε x, ε y, ε zBe respectively IMU coordinate system ox s, oy s, oz sThe constant value drift of axle gyro; a x, a yBe respectively IMU coordinate system ox s, oy sThe white noise error of axis accelerometer; δ K Gx, δ K Gy, δ K GzBe respectively IMU coordinate system ox s, oy s, oz sThe constant multiplier error of axle gyro; ω x, ω y, ω zBe respectively IMU coordinate system ox s, oy s, oz sThe white noise error of axle gyro;
The state-transition matrix of system is:
A = F 2 × 2 1 f 2 × 3 T ~ 2 × 2 O 2 × 6 F 3 × 2 2 F 3 × 3 3 O 3 × 2 T 3 × 6 O 8 × 2 O 8 × 3 O 8 × 2 O 8 × 6
F 2 × 2 1 = V N R n tan L 2 ω ie sin L + V E R n tan L - ( 2 ω ie sin L + 2 V E R n tan L ) 0
F 3 × 2 2 = 0 - 1 R m 1 R n 0 tan L R n 0
F 3 × 3 3 = 0 ω ie sin L + V E tan L R n - ( ω ie cos L + V E R n ) - ( ω ie sin L + V E tan L R n ) 0 - V N R m ω ie cos L + V E R n V N R m 0
f 2 × 3 = 0 - f U f N f U 0 f E
T ~ 2 × 2 = T 11 T 12 T 21 T 22
T 3 × 6 = - T 11 - T 12 - T 13 - T 11 ω x - T 12 ω y - T 13 ω z - T 21 - T 22 - T 23 - T 21 ω x - T 22 ω y - T 23 ω z - T 31 - T 32 - T 33 - T 31 ω x - T 32 ω y - T 33 ω z
V E, V NThe speed of representing east orientation, north orientation respectively; ω x, ω y, ω zThree input angular velocities representing gyro respectively; ω IeThe expression rotational-angular velocity of the earth; R m, R nRepresent earth meridian, fourth of the twelve Earthly Branches radius-of-curvature at the tenth of the twelve Earthly Branches respectively; L represents local latitude; f E, f N, f UBe expressed as respectively navigation coordinate system down east orientation, north orientation, day to specific force;
2) set up the measurement equation of Kalman filtering:
The measurement equation of describing the rotation strapdown inertial navitation system (SINS) with linear first-order differential equation is as follows:
Z(t)=H(t)X(t)+V(t)
Wherein: the measurement vector of etching system during Z (t) expression t; The measurement matrix of H (t) expression system; The measurement noise of V (t) expression system;
The system measurements matrix is:
H = 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
Amount is measured as the speed that obtains after the high-pass filtering:
Z = V ~ E V ~ N T .
CN2009100732422A 2009-11-20 2009-11-20 Rotating strapdown system on-site proving method based on digital high-passing filtering Expired - Fee Related CN101706287B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN2009100732422A CN101706287B (en) 2009-11-20 2009-11-20 Rotating strapdown system on-site proving method based on digital high-passing filtering

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN2009100732422A CN101706287B (en) 2009-11-20 2009-11-20 Rotating strapdown system on-site proving method based on digital high-passing filtering

Publications (2)

Publication Number Publication Date
CN101706287A true CN101706287A (en) 2010-05-12
CN101706287B CN101706287B (en) 2012-01-04

Family

ID=42376524

Family Applications (1)

Application Number Title Priority Date Filing Date
CN2009100732422A Expired - Fee Related CN101706287B (en) 2009-11-20 2009-11-20 Rotating strapdown system on-site proving method based on digital high-passing filtering

Country Status (1)

Country Link
CN (1) CN101706287B (en)

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101893445A (en) * 2010-07-09 2010-11-24 哈尔滨工程大学 Rapid initial alignment method for low-accuracy strapdown inertial navigation system under swinging condition
CN102175095A (en) * 2011-03-02 2011-09-07 浙江大学 Strap-down inertial navigation transfer alignment algorithm parallel implementation method
CN102435193A (en) * 2011-12-07 2012-05-02 浙江大学 High-precision initial alignment method of strapdown inertial navigation system
CN102538789A (en) * 2011-12-09 2012-07-04 北京理工大学 Rotating method of modulation type inertial navigation system with double-axis rotating continuously
CN102620735A (en) * 2012-04-17 2012-08-01 华中科技大学 Transposition method of double-shaft rotating type strapdown inertial navigation system for ship
CN102788596A (en) * 2012-08-16 2012-11-21 辽宁工程技术大学 Spot calibration method of rotary strap-down inertial navigation system with unknown carrier attitude
WO2013131471A1 (en) * 2012-03-06 2013-09-12 武汉大学 Quick calibration method for inertial measurement unit
CN103453917A (en) * 2013-09-04 2013-12-18 哈尔滨工程大学 Initial alignment and self-calibration method of double-shaft rotation type strapdown inertial navigation system
CN103604442A (en) * 2013-11-14 2014-02-26 哈尔滨工程大学 Observability analysis method applied to online calibration of strapdown inertial navitation system
CN103852085A (en) * 2014-03-26 2014-06-11 北京航空航天大学 Field calibration method of optical strapdown inertial navigation system based on least square fit
CN103852086A (en) * 2014-03-26 2014-06-11 北京航空航天大学 Field calibration method of optical fiber strapdown inertial navigation system based on Kalman filtering
CN104483064A (en) * 2014-12-29 2015-04-01 中国海洋石油总公司 In-situ calibration method for soft steel arm type mooring system stress monitoring device
CN104634364A (en) * 2015-01-29 2015-05-20 哈尔滨工程大学 Fiber-optic gyroscope scale factor self-calibration system based on step pulse modulation
CN106017507A (en) * 2016-05-13 2016-10-12 北京航空航天大学 Method for fast calibration of medium-and-low-precision optical fiber inertia units
CN110501027A (en) * 2019-09-16 2019-11-26 哈尔滨工程大学 A kind of optimal turn for dual-axis rotation MEMS-SINS stops time allocation method used therein

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR2725026B1 (en) * 1994-09-28 1997-01-10 Aerospatiale METHOD AND DEVICE FOR MINIMIZING IN AN INERTIAL MEASUREMENT SYSTEM THE ERROR DUE TO A MOVEMENT DISTURBING IN THE SPEED RESTITUTION
CN100554884C (en) * 2008-02-28 2009-10-28 北京航空航天大学 Method for standardization of optimum 8 positions of flexure gyroscope
CN101246022B (en) * 2008-03-21 2010-06-09 哈尔滨工程大学 Optic fiber gyroscope strapdown inertial navigation system two-position initial alignment method based on filtering
CN101377422B (en) * 2008-09-22 2010-09-08 北京航空航天大学 Method for calibrating optimum 24 positions of flexible gyroscope static drift error model
CN101514900B (en) * 2009-04-08 2011-01-26 哈尔滨工程大学 Method for initial alignment of a single-axis rotation strap-down inertial navigation system (SINS)
CN101514899B (en) * 2009-04-08 2010-12-01 哈尔滨工程大学 Optical fibre gyro strapdown inertial navigation system error inhibiting method based on single-shaft rotation
CN101571394A (en) * 2009-05-22 2009-11-04 哈尔滨工程大学 Method for determining initial attitude of fiber strapdown inertial navigation system based on rotating mechanism

Cited By (22)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101893445A (en) * 2010-07-09 2010-11-24 哈尔滨工程大学 Rapid initial alignment method for low-accuracy strapdown inertial navigation system under swinging condition
CN102175095A (en) * 2011-03-02 2011-09-07 浙江大学 Strap-down inertial navigation transfer alignment algorithm parallel implementation method
CN102175095B (en) * 2011-03-02 2013-06-19 浙江大学 Strap-down inertial navigation transfer alignment algorithm parallel implementation method
CN102435193A (en) * 2011-12-07 2012-05-02 浙江大学 High-precision initial alignment method of strapdown inertial navigation system
CN102538789A (en) * 2011-12-09 2012-07-04 北京理工大学 Rotating method of modulation type inertial navigation system with double-axis rotating continuously
CN102538789B (en) * 2011-12-09 2014-07-02 北京理工大学 Rotating method of modulation type inertial navigation system with double-axis rotating continuously
WO2013131471A1 (en) * 2012-03-06 2013-09-12 武汉大学 Quick calibration method for inertial measurement unit
CN102620735B (en) * 2012-04-17 2013-12-25 华中科技大学 Transposition method of double-shaft rotating type strapdown inertial navigation system for ship
CN102620735A (en) * 2012-04-17 2012-08-01 华中科技大学 Transposition method of double-shaft rotating type strapdown inertial navigation system for ship
CN102788596B (en) * 2012-08-16 2015-04-01 辽宁工程技术大学 Spot calibration method of rotary strap-down inertial navigation system with unknown carrier attitude
CN102788596A (en) * 2012-08-16 2012-11-21 辽宁工程技术大学 Spot calibration method of rotary strap-down inertial navigation system with unknown carrier attitude
CN103453917A (en) * 2013-09-04 2013-12-18 哈尔滨工程大学 Initial alignment and self-calibration method of double-shaft rotation type strapdown inertial navigation system
CN103604442A (en) * 2013-11-14 2014-02-26 哈尔滨工程大学 Observability analysis method applied to online calibration of strapdown inertial navitation system
CN103852085A (en) * 2014-03-26 2014-06-11 北京航空航天大学 Field calibration method of optical strapdown inertial navigation system based on least square fit
CN103852086A (en) * 2014-03-26 2014-06-11 北京航空航天大学 Field calibration method of optical fiber strapdown inertial navigation system based on Kalman filtering
CN103852085B (en) * 2014-03-26 2016-09-21 北京航空航天大学 A kind of fiber strapdown inertial navigation system system for field scaling method based on least square fitting
CN103852086B (en) * 2014-03-26 2016-11-23 北京航空航天大学 A kind of fiber strapdown inertial navigation system system for field scaling method based on Kalman filtering
CN104483064A (en) * 2014-12-29 2015-04-01 中国海洋石油总公司 In-situ calibration method for soft steel arm type mooring system stress monitoring device
CN104634364A (en) * 2015-01-29 2015-05-20 哈尔滨工程大学 Fiber-optic gyroscope scale factor self-calibration system based on step pulse modulation
CN106017507A (en) * 2016-05-13 2016-10-12 北京航空航天大学 Method for fast calibration of medium-and-low-precision optical fiber inertia units
CN106017507B (en) * 2016-05-13 2019-01-08 北京航空航天大学 A kind of used group quick calibrating method of the optical fiber of precision low used in
CN110501027A (en) * 2019-09-16 2019-11-26 哈尔滨工程大学 A kind of optimal turn for dual-axis rotation MEMS-SINS stops time allocation method used therein

Also Published As

Publication number Publication date
CN101706287B (en) 2012-01-04

Similar Documents

Publication Publication Date Title
CN101706287B (en) Rotating strapdown system on-site proving method based on digital high-passing filtering
CN101514899B (en) Optical fibre gyro strapdown inertial navigation system error inhibiting method based on single-shaft rotation
CN101514900B (en) Method for initial alignment of a single-axis rotation strap-down inertial navigation system (SINS)
CN103090867B (en) Error restraining method for fiber-optic gyroscope strapdown inertial navigation system rotating relative to geocentric inertial system
CN104655131B (en) Inertial navigation Initial Alignment Method based on ISTSSRCKF
CN103090870B (en) Spacecraft attitude measurement method based on MEMS (micro-electromechanical systems) sensor
CN103575299B (en) Utilize dual-axis rotation inertial navigation system alignment and the error correcting method of External Observation information
CN101713666B (en) Single-shaft rotation-stop scheme-based mooring and drift estimating method
CN101706284B (en) Method for increasing position precision of optical fiber gyro strap-down inertial navigation system used by ship
CN103090866B (en) Method for restraining speed errors of single-shaft rotation optical fiber gyro strapdown inertial navigation system
CN104374388B (en) Flight attitude determining method based on polarized light sensor
CN103900608B (en) A kind of low precision inertial alignment method based on quaternary number CKF
CN100547352C (en) The ground speed testing methods that is suitable for fiber optic gyro strapdown inertial navigation system
CN101629826A (en) Coarse alignment method for fiber optic gyro strapdown inertial navigation system based on single axis rotation
CN101571394A (en) Method for determining initial attitude of fiber strapdown inertial navigation system based on rotating mechanism
CN106153073B (en) A kind of nonlinear initial alignment method of full posture Strapdown Inertial Navigation System
CN103076025B (en) A kind of optical fibre gyro constant error scaling method based on two solver
CN105628025B (en) A kind of constant speed offset frequency/machine laser gyroscope shaking inertial navigation system air navigation aid
CN103697878B (en) A kind of single gyro list accelerometer rotation modulation north finding method
CN103076026B (en) A kind of method determining Doppler log range rate error in SINS
CN101881619A (en) Ship's inertial navigation and astronomical positioning method based on attitude measurement
CN106441357A (en) Damping network based single-axial rotation SINS axial gyroscopic drift correction method
CN103743413A (en) Installation error online estimation and north-seeking error compensation method for modulating north seeker under inclined state
CN106940193A (en) A kind of ship self adaptation based on Kalman filter waves scaling method
CN102116634A (en) Autonomous dimensionality reduction navigation method for deep sky object (DSO) landing detector

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20120104

Termination date: 20171120

CF01 Termination of patent right due to non-payment of annual fee