CN106996778B - Error parameter scaling method and device - Google Patents
Error parameter scaling method and device Download PDFInfo
- Publication number
- CN106996778B CN106996778B CN201710170915.0A CN201710170915A CN106996778B CN 106996778 B CN106996778 B CN 106996778B CN 201710170915 A CN201710170915 A CN 201710170915A CN 106996778 B CN106996778 B CN 106996778B
- Authority
- CN
- China
- Prior art keywords
- matrix
- error
- indicate
- vector
- indicates
- 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.)
- Active
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/10—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
- G01C21/12—Navigation; 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/16—Navigation; 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
Abstract
The invention discloses a kind of error parameter scaling method and devices.It include error parameter vector in the state equation and measurement equation this method comprises: establishing the state equation and measurement equation of aircraft guidance system, the error parameter vector is made of multiple error parameters；Differentiate the observability of each error parameter；When there are observable error parameter, using preset duration as filtering cycle, using Kalman filtering or adaptive-filtering, the observable error parameter is calibrated.The present invention realizes the purpose that real-time, on-orbit calibration goes out the error parameter of Guidance instrumentation.
Description
Technical field
The present invention relates to aircraft navigation control technology more particularly to a kind of error parameter scaling methods and device.
Background technique
Inertia type instrument is the heart of Strapdown Inertial Navigation System, and the size of error will directly affect the essence of entering the orbit of spacecraft
The size of degree and offset landings.Some spacecrafts for executing deep space exploration cannot reuse group after flying to certain altitude
Navigation is closed to be modified；Or important strike weapon, behind the native country that flies out, in order to improve its anti-interference ability, Neng Gouzheng
Combinations of satellites navigation feature can also be closed by really executing strike task.At this point, if being able to carry out an inertia system Guidance instrumentation
Error coefficient on-orbit calibration will be an important means of raising spacecraft navigation accuracy, and error coefficient mask work
Core be method for parameter estimation research.
Currently, the method comparative maturity of inertia system error parameter ground calibration, but inertia system practical application
When into aerial mission, due to being shaken by aircraft, the variation of space environment, the error parameter of ground calibration is unable to satisfy winged
Row device, spacecraft accurately determine the requirement of appearance positioning in space flight for a long time.
Summary of the invention
Technical problem solved by the present invention is compared with the prior art, a kind of error parameter scaling method and dress are provided
It sets, realizes the purpose that real-time, on-orbit calibration goes out the error parameter of Guidance instrumentation.
Above-mentioned purpose of the invention is achieved by the following technical programs:
In a first aspect, the present invention provides a kind of error parameter scaling methods, comprising:
The state equation and measurement equation of aircraft guidance system are established, includes in the state equation and measurement equation
Error parameter vector, the error parameter vector are made of multiple error parameters；
Differentiate the observability of each error parameter；
When there are observable error parameter, using preset duration as filtering cycle, using Kalman filtering or adaptively
Filtering, calibrates the observable error parameter.
Further, the error parameter includes the misaligned angle of the platform error, velocity error, location error, the scale of gyro
Factor error, gyroscopic drift error, the drift error of the scale factor error of accelerometer and accelerometer；It is each described
Error parameter is vector, and includes three durection components.
Further, the state equation are as follows:
In formula (1), A indicates state matrix；For the transition matrix of aircraft body system to launching inertial system；X indicates to miss
Poor parameter vector,φ indicates the misaligned angle of the platform error, δ
V indicates velocity error, and δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ Ka
Indicate the scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro
The white noise of measurement, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The shape
State matrix A are as follows:
In formula (2), The aircraft body system measured for gyro is relative to launching inertial system
Angular speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,Expression is asked
SolutionAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；G
Indicate that Newtonian gravitational constant, M indicate that earth quality, x, y, z indicate that coordinate of the aircraft under launching inertial system, r indicate flight
Distance of the device to launching inertial system origin；For the transition matrix of aircraft body system to launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi″For the attitude angle and inertial reference calculation posture of star sensor measurement
The difference at angle；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；ZrIt (t) is the aircraft position of GPS measurement
Set the difference with the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise to
Amount；The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor GPS's
Measure speed white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation.
Further, differentiate the observability of each error parameter, comprising:
Using the state matrix, transfer matrix, the calculation formula of the transfer matrix are calculated are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 list
Bit matrix；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2；
According to the transfer matrix and the measurement matrix, calculate to singular value decomposition matrix, it is described to singular value decomposition
The calculation formula of matrix are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,00 moment of expression state is to shape
The transfer matrix at 1 moment of state, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition fortune
It calculates；
Described singular value decomposition will be carried out to singular value decomposition matrix, to obtain the first singular value vector, the second singular value
Vector and singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i
∈[1,m]；U indicates the first singular value vector, U=[ui]=[u1u2...um]；V indicates the second singular value vector, V=[vi]=
[v1 v2 ... vm]；T indicates transposition operation；
Utilize the σi、ui、viAnd the guidance matrix of differences Z (t), calculate observability discriminant vector Y21×1, institute
State the calculation formula of observability discriminant vector Y are as follows:
In formula (9), T indicates transposition operation；
By the observability discriminant vector Y21×1Each element be compared with given threshold δ e；
As the observability discriminant vector Y21×1L, the element at the 1st column position of l ∈ [1,21] row is greater than setting
When threshold value δ e, l in decision errors parameter vector X, the error parameter Observable at the 1st column position of l ∈ [1,21] row.
Second aspect, the present invention also provides a kind of error parameter caliberating device, which includes:
Module is constructed, for establishing the state equation and measurement equation of aircraft guidance system, the state equation and amount
Surveying includes error parameter vector in equation, and the error parameter vector is made of multiple error parameters；
Discrimination module, for differentiating the observability of each error parameter；
Demarcating module, for using preset duration as filtering cycle, utilizing Kalman when there are observable error parameter
Filtering or adaptive-filtering, calibrate the observable error parameter.
Further, the error parameter includes the misaligned angle of the platform error, velocity error, location error, the scale of gyro
Factor error, gyroscopic drift error, the drift error of the scale factor error of accelerometer and accelerometer；It is each described
Error parameter is vector, and includes three durection components.
Further, the state equation are as follows:
In formula (1), A indicates state matrix；Cb iFor the transition matrix of aircraft body system to launching inertial system；X indicates to miss
Poor parameter vector,φ indicates the misaligned angle of the platform error, δ
V indicates velocity error, and δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ Ka
Indicate the scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro
The white noise of measurement, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The shape
State matrix A are as follows:
In formula (2), The aircraft body system measured for gyro is relative to launching inertial system
Angular speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,Expression is asked
SolutionAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；
G expression Newtonian gravitational constant, M expression earth quality, x,y, coordinate of the z expression aircraft under launching inertial system, r expression flight
Distance of the device to launching inertial system origin；For the transition matrix of aircraft body system to launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi″For the attitude angle and inertial reference calculation posture of star sensor measurement
The difference at angle；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；ZrIt (t) is the aircraft position of GPS measurement
Set the difference with the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise to
Amount；The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor the amount of GPS
Degree of testing the speed white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation.
Further, the discrimination module includes:
First computing unit calculates transfer matrix, the calculation formula of the transfer matrix for utilizing the state matrix
Are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 list
Bit matrix；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2；
Second computing unit, for calculating to singular value decomposition matrix according to the transfer matrix and the measurement matrix,
The calculation formula to singular value decomposition matrix are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,00 moment of expression state is to shape
The transfer matrix at 1 moment of state, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition fortune
It calculates；
Decomposition unit, for described singular value decomposition will to be carried out to singular value decomposition matrix, with obtain the first singular value to
Amount, the second singular value vector and singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i
∈[1,m]；U indicates the first singular value vector, U=[ui]=[u1 u2 ... um]；V indicates the second singular value vector, V=[vi]
=[v1 v2 ... vm]；T indicates transposition operation；
Third computing unit, for utilizing the σi、ui、viAnd the guidance matrix of differences Z (t), calculate Observable
Property discriminant vector Y21×1, the calculation formula of the observability discriminant vector Y are as follows:
In formula (9), T indicates transposition operation；
Comparing unit is used for the observability discriminant vector Y21×1Each element compared with given threshold δ e
Compared with；
Judging unit, for working as the observability discriminant vector Y21×1L, at the 1st column position of l ∈ [1,21] row
When element is greater than given threshold δ e, l in decision errors parameter vector X, the error at the 1st column position of l ∈ [1,21] row
Parameter Observable.
Compared with prior art, the present invention has the following advantages:
(1), the state equation and measurement equation of the invention by establishing aircraft guidance system；Differentiate each error parameter
Observability；When there are observable error parameter, using Kalman filtering or adaptive-filtering, calibrate described considerable
The error parameter of survey can overcome the shortcomings of existing guidance instrument error isolation technics " world is inconsistent ", can be real-time, in-orbit
The error coefficient of Guidance instrumentation is calibrated, algorithm is simple, convenient for engineering.
(2), the present invention can reduce inertia device Support expense, due to the ability with on-orbit calibration, so as to
The number for enough reducing inertia system ground periodic calibrating, saves human and material resources.
(3), the present invention can be improved aircraft guidance precision, reduce the dependence to satellite navigation, extending space aircraft
The range and ability of execution task.
Detailed description of the invention
Fig. 1 is the flow chart of one of embodiment of the present invention one error parameter scaling method；
Fig. 2 is the structure chart of one of embodiment of the present invention two error parameter caliberating device.
Specific embodiment
Invention is further described in detail with reference to the accompanying drawings and examples.It is understood that described herein
Specific embodiment be used only for explaining the present invention rather than limiting the invention.It also should be noted that for the ease of
It describes, only the parts related to the present invention are shown rather than entire infrastructure in attached drawing.
Embodiment one
Fig. 1 is the flow chart of one of embodiment of the present invention one error parameter scaling method, and the present embodiment is applicable to
The case where on-orbit calibration is carried out to the error parameter of Guidance instrumentation is needed, this method can be held by error parameter caliberating device
Row, wherein the device can be by software and or hardware realization.With reference to Fig. 1, error parameter scaling method tool provided in this embodiment
Body may include steps of:
S110, the state equation and measurement equation for establishing aircraft guidance system, in the state equation and measurement equation
It include error parameter vector, the error parameter vector is made of multiple error parameters.
Specifically, the error parameter include the misaligned angle of the platform error, velocity error, location error, the scale of gyro because
Number error, gyroscopic drift error, the drift error of the scale factor error of accelerometer and accelerometer；Each mistake
Poor parameter is vector, and includes three durection components.Each error parameter is vector, and includes three directions
Component, i.e., each error parameter is vector, includes three direction vectors in respective coordinate system.For example, gyroscopic drift error
Include three durection components in its coordinate system: x-axis is to component, y-axis to component and z-axis to component.In another example accelerometer
Scale factor error in its coordinate system include three durection components: x-axis is to component, y-axis to component and z-axis to component.
Specifically, the state equation are as follows:
In formula (1), A indicates state matrix；For the transition matrix of aircraft body system to launching inertial system；X indicates to miss
Poor parameter vector,φ indicates the misaligned angle of the platform error, δ
V indicates velocity error, and δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ Ka
Indicate the scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro
The white noise of measurement, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The shape
State matrix A are as follows:
In formula (2), The aircraft body system measured for gyro is relative to launching inertial system
Angular speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,It indicates
It solvesAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；
G indicates that Newtonian gravitational constant, M indicate that earth quality, x, y, z indicate that coordinate of the aircraft under launching inertial system, r indicate flight
Distance of the device to launching inertial system origin；For the transition matrix of aircraft body system to launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi" it is attitude angle and inertial reference calculation posture that star sensor measures
The difference at angle；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；ZrIt (t) is the aircraft position of GPS measurement
Set the difference with the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise to
Amount；The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor the amount of GPS
Degree of testing the speed white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation.
S120, the observability for differentiating each error parameter.
Optionally, differentiate the observability of each error parameter, comprising:
Using the state matrix, transfer matrix, the calculation formula of the transfer matrix are calculated are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 list
Bit matrix；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2.
Specifically, T1Indicate preset duration, the preset duration is usually the integral multiple of 200ms, is up to 1s.For example,
If preset duration T1=200ms, AkIndicate the state matrix of k-th of 200ms period.
According to the transfer matrix and the measurement matrix, calculate to singular value decomposition matrix, it is described to singular value decomposition
The calculation formula of matrix are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,00 moment of expression state is to shape
The transfer matrix at 1 moment of state, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition fortune
It calculates.
Described singular value decomposition will be carried out to singular value decomposition matrix, to obtain the first singular value vector, the second singular value
Vector and singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i
∈[1,m]；U indicates the first singular value vector, U=[ui]=[u1 u2 ... um]；V indicates the second singular value vector, V=[vi]
=[v1 v2 ... vm]；T indicates transposition operation.
Utilize the σi、ui、viAnd the guidance matrix of differences Z (t), calculate observability discriminant vector Y21×1, institute
State the calculation formula of observability discriminant vector Y are as follows:
In formula (9), T indicates transposition operation.
By the observability discriminant vector Y21×1Each element be compared with given threshold δ e.
Specifically, due to observability discriminant vector Y21×1In each element range between 0~1, it is thus, described
Given threshold δ e is chosen as 0.1.In the present embodiment, the given threshold δ e=0.1, i.e., by the observability discriminant vector
Y21×1Each element be compared with given threshold 0.1.
As the observability discriminant vector Y21×1L, the element at the 1st column position of l ∈ [1,21] row is greater than setting
When threshold value δ e, l in decision errors parameter vector X, the error parameter Observable at the 1st column position of l ∈ [1,21] row.
Specifically, in the present embodiment,X indicates to miss
Poor parameter vector, φ indicate the misaligned angle of the platform error, and δ V indicates velocity error, and δ r indicates location error, δ KgIndicate the mark of gyro
Spend factor error, b1Indicate gyroscopic drift error, δ KaIndicate the scale factor error of accelerometer,Indicate accelerometer
Drift error；And each error parameter includes three direction vectors in respective coordinate system.For example, when the observability differentiates
Vector Y21×1The 1st column position of the 10th row at element when being greater than given threshold δ e=0.1, the in decision errors parameter vector X
The error parameter Observable at the 1st column position of 10 row, i.e. the scale factor error δ K of gyrogObservable.In another example working as institute
State observability discriminant vector Y21×1The 1st column position of the 19th row at element be greater than given threshold δ e=0.1 when, decision errors
The error parameter Observable in parameter vector X at the 1st column position of the 19th row, the i.e. drift error of accelerometerIt is considerable
It surveys.
S130, when there are observable error parameter, using preset duration as filtering cycle, using Kalman filtering or from
Adaptive filtering calibrates the observable error parameter.
Specifically, in the present embodiment, when there are observable error parameter, with the preset duration in step S120
T1For filtering cycle, the filtering cycle can be taken as the integral multiple of 200ms, be up to 1s.Using Kalman filtering or adaptively
Filtering, calibrates the observable error parameter.
The technical solution of the present embodiment is by establishing the state equation and measurement equation of aircraft guidance system；Differentiate each
The observability of error parameter；When there are observable error parameter, using Kalman filtering or adaptive-filtering, calibrate
The observable error parameter can overcome the shortcomings of existing guidance instrument error isolation technics " world is inconsistent ", can
In real time, on-orbit calibration goes out the error coefficient of Guidance instrumentation, and algorithm is simple, convenient for engineering；It can reduce inertia device maintenance to protect
Barrier expense, since the ability with on-orbit calibration saves manpower so as to reduce the number of inertia system ground periodic calibrating
Material resources；It can be improved aircraft guidance precision, reduce the dependence to satellite navigation, extending space aircraft executes the range of task
And ability.
Embodiment two
Fig. 2 is the structure chart of one of embodiment of the present invention two error parameter caliberating device, and the present embodiment is applicable to
Need the case where on-orbit calibration is carried out to the error parameter of Guidance instrumentation.With reference to Fig. 2, error parameter calibration provided in this embodiment
Device specifically can be such that
Construct module 210, for establishing the state equation and measurement equation of aircraft guidance system, the state equation and
It include error parameter vector in measurement equation, the error parameter vector is made of multiple error parameters；
Discrimination module 220, for differentiating the observability of each error parameter；
Demarcating module 230, for using preset duration as filtering cycle, utilizing card when there are observable error parameter
Kalman Filtering or adaptive-filtering calibrate the observable error parameter.
Optionally, the error parameter include the misaligned angle of the platform error, velocity error, location error, the scale of gyro because
Number error, gyroscopic drift error, the drift error of the scale factor error of accelerometer and accelerometer；Each mistake
Poor parameter is vector, and includes three durection components.
Optionally, the state equation are as follows:
In formula (1), A indicates state matrix；For the transition matrix of aircraft body system to launching inertial system；X indicates to miss
Poor parameter vector,φ indicates the misaligned angle of the platform error, δ
V indicates velocity error, and δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ Ka
Indicate the scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro
The white noise of measurement, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The shape
State matrix A are as follows:
In formula (2), The aircraft body system measured for gyro is relative to launching inertial system
Angular speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,Table
Show solutionAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；
G indicates that Newtonian gravitational constant, M indicate that earth quality, x, y, z indicate that coordinate of the aircraft under launching inertial system, r indicate flight
Distance of the device to launching inertial system origin；For the transition matrix of aircraft body system to launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi" it is attitude angle and inertial reference calculation posture that star sensor measures
The difference at angle；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；ZrIt (t) is the aircraft position of GPS measurement
Set the difference with the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise to
Amount；The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor the amount of GPS
Degree of testing the speed white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation.
Optionally, the discrimination module includes:
First computing unit calculates transfer matrix, the calculation formula of the transfer matrix for utilizing the state matrix
Are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 list
Bit matrix；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2；
Second computing unit, for calculating to singular value decomposition matrix according to the transfer matrix and the measurement matrix,
The calculation formula to singular value decomposition matrix are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,00 moment of expression state is to shape
The transfer matrix at 1 moment of state, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition fortune
It calculates；
Decomposition unit, for described singular value decomposition will to be carried out to singular value decomposition matrix, with obtain the first singular value to
Amount, the second singular value vector and singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i
∈[1,m]；U indicates the first singular value vector, U=[ui]=[u1 u2 ... um]；V indicates the second singular value vector, V=[vi]
=[v1 v2 ... vm]；T indicates transposition operation；
Third computing unit, for utilizing the σi、ui、viAnd the guidance matrix of differences Z (t), calculate Observable
Property discriminant vector Y21×1, the calculation formula of the observability discriminant vector Y are as follows:
In formula (9), T indicates transposition operation；
Comparing unit is used for the observability discriminant vector Y21×1Each element compared with given threshold δ e
Compared with；
Judging unit, for working as the observability discriminant vector Y21×1L, at the 1st column position of l ∈ [1,21] row
When element is greater than given threshold δ e, l in decision errors parameter vector X, the error at the 1st column position of l ∈ [1,21] row
Parameter Observable.
Error parameter caliberating device provided in this embodiment is demarcated with error parameter provided by any embodiment of the invention
Method belongs to same inventive concept, and error parameter scaling method provided by any embodiment of the invention can be performed, have execution
The corresponding functional module of error parameter scaling method and beneficial effect.The not technical detail of detailed description in the present embodiment, can
The error parameter scaling method provided referring to any embodiment of that present invention.
Note that the above is only a better embodiment of the present invention and the applied technical principle.It will be appreciated by those skilled in the art that
The invention is not limited to the specific embodiments described herein, be able to carry out for a person skilled in the art it is various it is apparent variation,
It readjusts and substitutes without departing from protection scope of the present invention.Therefore, although being carried out by above embodiments to the present invention
It is described in further detail, but the present invention is not limited to the above embodiments only, without departing from the inventive concept, also
It may include more other equivalent embodiments, and the scope of the invention is determined by the scope of the appended claims.
Claims (3)
1. a kind of error parameter scaling method characterized by comprising
The state equation and measurement equation of aircraft guidance system are established, includes error in the state equation and measurement equation
Parameter vector, the error parameter vector are made of multiple error parameters；
Differentiate the observability of each error parameter；
When there are observable error parameter, using preset duration as filtering cycle, using Kalman filtering or adaptive-filtering,
Calibrate the observable error parameter, in which:
The state equation are as follows:
In formula (1), A indicates state matrix；For the transition matrix of aircraft body system to launching inertial system；X indicates error ginseng
Number vector,φ indicates the misaligned angle of the platform error, δ V table
Show velocity error, δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ KaIt indicates
The scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro to measure
White noise, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The state square
Battle array A are as follows:
In formula (2), Angle of the aircraft body system measured for gyro relative to launching inertial system
Speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,Table
Show solutionAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；G indicates that Newtonian gravitational constant, M indicate earth quality, x, y, z table
Show that coordinate of the aircraft under launching inertial system, r indicate aircraft to the distance of launching inertial system origin；For aircraft sheet
Transition matrix of the system to launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi”For star sensor measurement attitude angle and inertial reference calculation attitude angle it
Difference；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；Zr(t) for GPS measurement position of aircraft and
The difference of the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise vector；
The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor the measurement speed of GPS
Spend white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation；
Differentiate the observability of each error parameter method particularly includes:
Using the state matrix, transfer matrix, the calculation formula of the transfer matrix are calculated are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 unit square
Battle array；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2；
According to the transfer matrix and the measurement matrix, calculate to singular value decomposition matrix, it is described to singular value decomposition matrix
Calculation formula are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,0When 0 moment of expression state is to state 1
The transfer matrix at quarter, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition operation；
Described will carry out singular value decomposition to singular value decomposition matrix, with obtain the first singular value vector, the second singular value vector,
And singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i∈[1,
m]；U indicates the first singular value vector, U=[ui]=[u1 u2 ... um]；V indicates the second singular value vector, V=[vi]=[v1
v2 ... vm]；T indicates transposition operation；
Utilize the σi、ui、viAnd the guidance matrix of differences Z (t), calculate observability discriminant vector Y21×1, it is described can
The calculation formula of observation discriminant vector Y are as follows:
In formula (9), T indicates transposition operation；
By the observability discriminant vector Y21×1Each element be compared with given threshold δ e；
As the observability discriminant vector Y21×1L, the element at the 1st column position of l ∈ [1,21] row is greater than given threshold δ
When e, l in decision errors parameter vector X, the error parameter Observable at the 1st column position of l ∈ [1,21] row；
The error parameter includes the misaligned angle of the platform error, velocity error, location error, the scale factor error of gyro, gyro
The drift error of drift error, the scale factor error of accelerometer and accelerometer；Each error parameter is arrow
Amount, and include three durection components.
2. a kind of error parameter caliberating device characterized by comprising
Module is constructed, for establishing the state equation and measurement equation of aircraft guidance system, the state equation and measurement side
Cheng Zhongjun includes error parameter vector, and the error parameter vector is made of multiple error parameters；
Discrimination module, for differentiating the observability of each error parameter；
The discrimination module includes:
First computing unit calculates transfer matrix, the calculation formula of transfer matrix for utilizing state matrix are as follows:
In formula (6), Φk,k-1Transfer matrix of the expression state k-1 moment to the state k moment, I21×21Indicate 21 × 21 unit square
Battle array；T1Indicate preset duration；AkIndicate the state matrix of k-th of preset duration；M is positive integer and m >=2；
Second computing unit, for calculating to singular value decomposition matrix, to singular value according to the transfer matrix and measurement matrix
The calculation formula of split-matrix are as follows:
In formula (7), Q is indicated to singular value decomposition matrix, H9×21Indicate measurement matrix, Φ1,0When 0 moment of expression state is to state 1
The transfer matrix at quarter, Φm-1,m-2The transfer matrix at expression state m-2 moment to state m-1 moment, T indicate transposition operation；
Decomposition unit, for described singular value decomposition will to be carried out to singular value decomposition matrix, to obtain the first singular value vector, the
Two singular value vectors and singular value matrix, the formula of the singular value decomposition are as follows:
Q=U ∑ VT (8)
In formula (8), Q is indicated to singular value decomposition matrix, and ∑ indicates that singular value matrix, the elements in a main diagonal of ∑ are σi,i∈[1,
m]；U indicates the first singular value vector, U=[ui]=[u1 u2 ... um]；V indicates the second singular value vector, V=[vi]=[v1
v2 ... vm]；T indicates transposition operation；
Third computing unit, for utilizing σi、ui、viAnd guidance matrix of differences Z (t), calculate observability discriminant vector
Y21×1, the calculation formula of observability discriminant vector Y are as follows:
In formula (9), T indicates transposition operation；
Comparing unit is used for the observability discriminant vector Y21×1Each element be compared with given threshold δ e；
Judging unit, for working as the observability discriminant vector Y21×1L, the element at the 1st column position of l ∈ [1,21] row
When greater than given threshold δ e, l in decision errors parameter vector X, the error parameter at the 1st column position of l ∈ [1,21] row
Observable；
Demarcating module, for using preset duration as filtering cycle, utilizing Kalman filtering when there are observable error parameter
Or adaptive-filtering, calibrate the observable error parameter, in which:
The state equation are as follows:
In formula (1), A indicates state matrix；For the transition matrix of aircraft body system to launching inertial system；X indicates error ginseng
Number vector,φ indicates the misaligned angle of the platform error, δ V table
Show velocity error, δ r indicates location error, δ KgIndicate the scale factor error of gyro, b1Indicate gyroscopic drift error, δ KaIt indicates
The scale factor error of accelerometer,Indicate the drift error of accelerometer, T indicates transposition operation；ηgIndicate gyro to measure
White noise, εaFor the white noise of accelerometer measures；For the first derivative vector of error parameter vector X；The state square
Battle array A are as follows:
In formula (2), Angle of the aircraft body system measured for gyro relative to launching inertial system
Speed aircraft body system projection,Indicate withFor the matrix of the elements in a main diagonal；fbThe specific force measured for accelerometer aircraft body system projection,Table
Show solutionAntisymmetric matrix；I3×3Indicate 3 × 3 unit matrix；G indicates that Newtonian gravitational constant, M indicate that earth quality, x, y, z indicate
Coordinate of the aircraft under launching inertial system, r indicate aircraft to the distance of launching inertial system origin；For aircraft body system
To the transition matrix of launching inertial system；
The measurement equation are as follows:
In formula (3), Z (t) indicates guidance matrix of differences；φi”For star sensor measurement attitude angle and inertial reference calculation attitude angle it
Difference；ZvIt (t) is the difference of the speed of the aircraft speed and inertial reference calculation of GPS measurement；Zr(t) for GPS measurement position of aircraft and
The difference of the position of inertial reference calculation；H9×21Indicate measurement matrix；X indicates the error parameter vector；V9×1Indicate white noise vector；
The measurement matrix H9×21Are as follows:
In formula (4), I3×3Indicate 3 × 3 unit matrix；The white noise vector V9×1Are as follows:
In formula (5),WithWhite noise, δ M are measured for the posture of star sensorx、δMyWith δ MzFor the measurement speed of GPS
Spend white noise, δ xG、δxGWith δ zGFor the adjustment location white noise of GPS；T indicates transposition operation.
3. the apparatus of claim 2, which is characterized in that the error parameter includes the misaligned angle of the platform error, speed
Error, location error, the scale factor error of gyro, gyroscopic drift error, the scale factor error of accelerometer and acceleration
Spend the drift error of meter；Each error parameter is vector, and includes three durection components.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710170915.0A CN106996778B (en) | 2017-03-21 | 2017-03-21 | Error parameter scaling method and device |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710170915.0A CN106996778B (en) | 2017-03-21 | 2017-03-21 | Error parameter scaling method and device |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106996778A CN106996778A (en) | 2017-08-01 |
| CN106996778B true CN106996778B (en) | 2019-11-29 |
Family
ID=59431707
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710170915.0A Active CN106996778B (en) | 2017-03-21 | 2017-03-21 | Error parameter scaling method and device |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN106996778B (en) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110608741A (en) * | 2019-09-25 | 2019-12-24 | 北京安达维尔科技股份有限公司 | Method for improving attitude calculation precision of aircraft attitude and heading reference system |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1995916A (en) * | 2006-12-27 | 2007-07-11 | 北京航空航天大学 | Integrated navigation method based on star sensor calibration |
| CN102297695A (en) * | 2010-06-22 | 2011-12-28 | 中国船舶重工集团公司第七○七研究所 | Kalman filtering method in deep integrated navigation system |
| CN102879011A (en) * | 2012-09-21 | 2013-01-16 | 北京控制工程研究所 | Lunar inertial navigation alignment method assisted by star sensor |
| CN103245359A (en) * | 2013-04-23 | 2013-08-14 | 南京航空航天大学 | Method for calibrating fixed errors of inertial sensor in inertial navigation system in real time |
| CN103557871A (en) * | 2013-10-22 | 2014-02-05 | 北京航空航天大学 | Strap-down inertial navigation air initial alignment method for floating aircraft |
| CN104764467A (en) * | 2015-04-08 | 2015-07-08 | 南京航空航天大学 | Online adaptive calibration method for inertial sensor errors of aerospace vehicle |
| CN105352529A (en) * | 2015-11-16 | 2016-02-24 | 南京航空航天大学 | Multisource-integrated-navigation-system distributed inertia node total-error on-line calibration method |
| CN105867399A (en) * | 2016-04-18 | 2016-08-17 | 北京航天自动控制研究所 | Method for determining multi-state tracking guidance parameters |
| CN106052716A (en) * | 2016-05-25 | 2016-10-26 | 南京航空航天大学 | Method for calibrating gyro errors online based on star light information assistance in inertial system |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7646336B2 (en) * | 2006-03-24 | 2010-01-12 | Containertrac, Inc. | Automated asset positioning for location and inventory tracking using multiple positioning techniques |
| US9417091B2 (en) * | 2013-05-13 | 2016-08-16 | The Johns Hopkins University | System and method for determining and correcting field sensors errors |
-
2017
- 2017-03-21 CN CN201710170915.0A patent/CN106996778B/en active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1995916A (en) * | 2006-12-27 | 2007-07-11 | 北京航空航天大学 | Integrated navigation method based on star sensor calibration |
| CN102297695A (en) * | 2010-06-22 | 2011-12-28 | 中国船舶重工集团公司第七○七研究所 | Kalman filtering method in deep integrated navigation system |
| CN102879011A (en) * | 2012-09-21 | 2013-01-16 | 北京控制工程研究所 | Lunar inertial navigation alignment method assisted by star sensor |
| CN103245359A (en) * | 2013-04-23 | 2013-08-14 | 南京航空航天大学 | Method for calibrating fixed errors of inertial sensor in inertial navigation system in real time |
| CN103557871A (en) * | 2013-10-22 | 2014-02-05 | 北京航空航天大学 | Strap-down inertial navigation air initial alignment method for floating aircraft |
| CN104764467A (en) * | 2015-04-08 | 2015-07-08 | 南京航空航天大学 | Online adaptive calibration method for inertial sensor errors of aerospace vehicle |
| CN105352529A (en) * | 2015-11-16 | 2016-02-24 | 南京航空航天大学 | Multisource-integrated-navigation-system distributed inertia node total-error on-line calibration method |
| CN105867399A (en) * | 2016-04-18 | 2016-08-17 | 北京航天自动控制研究所 | Method for determining multi-state tracking guidance parameters |
| CN106052716A (en) * | 2016-05-25 | 2016-10-26 | 南京航空航天大学 | Method for calibrating gyro errors online based on star light information assistance in inertial system |
Also Published As
| Publication number | Publication date |
|---|---|
| CN106996778A (en) | 2017-08-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106989761B (en) | A kind of spacecraft Guidance instrumentation on-orbit calibration method based on adaptive-filtering | |
| CN106289246B (en) | A kind of flexible link arm measure method based on position and orientation measurement system | |
| Zheng et al. | An eight-position self-calibration method for a dual-axis rotational inertial navigation system | |
| CN103913181B (en) | A kind of airborne distributed POS Transfer Alignments based on parameter identification | |
| CN110887507B (en) | Method for quickly estimating all zero offsets of inertial measurement units | |
| CN107101649B (en) | A kind of in-orbit error separating method of spacecraft Guidance instrumentation | |
| CN109724624B (en) | Airborne self-adaptive transfer alignment method suitable for wing deflection deformation | |
| CN108548542B (en) | Near-earth orbit determination method based on atmospheric resistance acceleration measurement | |
| Yao et al. | Transverse Navigation under the Ellipsoidal Earth Model and its Performance in both Polar and Non-polar areas | |
| Lu et al. | An all-parameter system-level calibration for stellar-inertial navigation system on ground | |
| CN107728182A (en) | Flexible more base line measurement method and apparatus based on camera auxiliary | |
| Rad et al. | Optimal attitude and position determination by integration of INS, star tracker, and horizon sensor | |
| CN107764261B (en) | Simulation data generation method and system for distributed POS (point of sale) transfer alignment | |
| CN104698486A (en) | Real-time navigation method of data processing computer system for distributed POS | |
| Gebre-Egziabher et al. | MAV attitude determination by vector matching | |
| CN109931955A (en) | Strapdown inertial navigation system Initial Alignment Method based on the filtering of state correlation Lie group | |
| CN110296719B (en) | On-orbit calibration method | |
| CN104833375B (en) | A kind of IMU Two position methods by star sensor | |
| CN109708663B (en) | Star sensor online calibration method based on aerospace plane SINS assistance | |
| Li et al. | Fast fine initial self-alignment of INS in erecting process on stationary base | |
| CN116105730A (en) | Angle measurement-only optical combination navigation method based on cooperative target satellite very short arc observation | |
| CN105606093B (en) | Inertial navigation method and device based on gravity real-Time Compensation | |
| CN109581523B (en) | Method and system for calibrating accelerometer by satellite tracking satellite device | |
| CN106996778B (en) | Error parameter scaling method and device | |
| Yuan et al. | Dynamic initial alignment of the MEMS-based low-cost SINS for AUV based on unscented Kalman filter |
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 |