CN107942090A - A kind of spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information - Google Patents

A kind of spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information Download PDF

Info

Publication number
CN107942090A
CN107942090A CN201711456224.3A CN201711456224A CN107942090A CN 107942090 A CN107942090 A CN 107942090A CN 201711456224 A CN201711456224 A CN 201711456224A CN 107942090 A CN107942090 A CN 107942090A
Authority
CN
China
Prior art keywords
mrow
msub
mtd
msup
spacecraft
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
CN201711456224.3A
Other languages
Chinese (zh)
Other versions
CN107942090B (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.)
Beihang University
Original Assignee
Beihang University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beihang University filed Critical Beihang University
Priority to CN201711456224.3A priority Critical patent/CN107942090B/en
Publication of CN107942090A publication Critical patent/CN107942090A/en
Application granted granted Critical
Publication of CN107942090B publication Critical patent/CN107942090B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P3/00—Measuring linear or angular speed; Measuring differences of linear or angular speeds
    • G01P3/36—Devices characterised by the use of optical means, e.g. using infrared, visible, or ultraviolet light

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Electromagnetism (AREA)
  • Power Engineering (AREA)
  • General Physics & Mathematics (AREA)
  • Navigation (AREA)
  • Image Analysis (AREA)
  • Feedback Control In General (AREA)

Abstract

The present invention relates to a kind of spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information.Selection is used as state model using spacecraft angular speed and Singer acceleration probabilistic model of the angular acceleration as quantity of state, star chart is obscured to multiframe and carries out image procossing, Optic flow information is extracted as measurement, measurement model is established according to the relation between asterism light stream and spacecraft angular speed, the angular speed of three axis of spacecraft can be estimated finally by Kalman filtering.The present invention is directed to the gyro free Attitude estimation problem under spacecraft high dynamic condition, in spacecraft rotation at a high speed, star sensor image blur, realizes and is directly estimated using fuzzy star chart not against the attitude angular velocity of importance in star map recognition process.

Description

A kind of spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information
Technical field
The invention belongs to spacecraft attitude to determine field, is related to a kind of spacecraft angle speed based on fuzzy star chart Optic flow information Spend method of estimation.
Background technology
For spacecraft when performing aerial mission, the attitude motion ginseng such as its direction, rotation angle and angular speed Number requires accurate.The precision that attitude motion of spacecraft parameter determines directly determines the precision of gesture stability, and finally determines The function of whole space mission is realized.
Existing Attitude rate estimator method is typically to rely on gyro output spacecraft attitude and angular velocity information, but its cost It is higher, and there is the problems such as drift and noise, therefore many gyro free Attitude rate estimator methods were proposed in recent years to solve this A little problems and for preventing gyro failure.Gyro free Attitude rate estimator method is to obtain measurement based on Sensitive Apparatus mostly, Such as star sensor, magnetometer etc..When carrying out survey of deep space when task, it is necessary to realize that the spacecraft of accurate three axis control is usual The star sensor of installation can be relied on to carry out posture to determine.The existing gyro free Attitude rate estimator method based on star sensor is to work as When spacecraft low-speed motion or shorter time for exposure, star sensor can obtain more clearly image, by pre- to star chart Processing, solves the attitude parameter of spacecraft after barycenter extraction and match cognization.It is longer when spacecraft high-speed motion or time for exposure When, it can make to be moved to one other pixel point, that is, an asterism pair from a pixel during star sensor progress asterism imaging Multiple pixels are answered, track is in ribbon, and acquisition is fuzzy star chart.Asterism centroid calculation precision can reduce at this time, sensitive to star The application of device causes limitation.
For obtaining the situation that image is fuzzy star chart, existing method be by its deblurring, after image recovery method, Barycenter extraction is carried out to picture rich in detail, posture solution is carried out after importance in star map recognition and matching, but directly obtained using fuzzy star chart Measurement information carry out attitude angular velocity estimation method almost without.
The content of the invention
The technical problem to be solved in the present invention is：When the star chart of acquisition is fuzzy star chart, without deblurring process, directly Connect and extract Optic flow information using fuzzy star chart, all inapplicable in Region Matching method and the optical flow computation method of approximate asterism displacement In the case of, a kind of computational methods for longer or non-linear asterism blurring trajectorie of design, image procossing is carried out to every frame star chart Extraction Optic flow information is estimated as measurement, and then to spacecraft angular speed afterwards.
The present invention proposes a kind of spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information, passes through selection Using spacecraft angular speed and angular acceleration as the Singer acceleration probabilistic model of quantity of state as state model, to m frames, m ＞ 2, obscures star chart and carries out image procossing, and using the Optic flow information of extraction as measurement, according to asterism light stream and three-axis stabilization Relation between spacecraft angular speed establishes measurement model, finally by the angle speed of three axis of Kalman Filter Estimation spacecraft Degree.Specifically include following steps：
1. the foundation of System State Model；
First, angular speed and angular acceleration of the quantity of state for spacecraft are selected, i.e.,：
In formula, ω1(k), ω2(k) and ω3(k) when being kth frame (1≤k≤m) star chart respectively spacecraft in camera coordinates Around x-axis under system, the angular speed of y-axis and z-axis,WithRespectively kth frame star chart when spacecraft in three axis Angular acceleration.State model selects Singer acceleration probabilistic models, and state equation is as follows：
X (k)=Φ (T) X (k-1)+W (k-1) (2)
In formula, T is filtering cycle, Λ=diag { 1/ τ1,1/τ2,1/τ3, wherein τi(i=1,2,3), represents spacecraft The G-time constant of the corresponding x-axis of body, y-axis and z-axis, design parameter value are set in Kalman filtering process, and I is unit Matrix；W (k-1) is estimated system noise sequence suffered by state X (k), the variance matrix representation of W (k-1) for the kth frame star chart moment For Q (k-1), its each element is expressed as：
Q11(k-1)=Λ-2Σ2(4e-ΛT-3I-e-2ΛT+2ΛT) (4)
Q21(k-1)=Q12(k-1)=Λ-1Σ2(e-2ΛT+I-2e-ΛT) (5)
Q22(k-1)=Σ2(I-e-2ΛT) (6)
In formula, Σ=diag { σ1,σ2,σ3It is used for determining the population variance of noise in Singer acceleration probabilistic models, It is expressed as：
Wherein, WithAngular acceleration element can be regarded asPdf model adjusting parameter, It is expressed as：Angular accelerationEqual to maximumAnd minimum valueProbability beProbability equal to 0 ForIn sectionOn be uniformly distributed, probability isParameter WithScope It is as follows：Occurrence is selected in Kalman filtering process.
2. measurement obtains and the foundation of system measurements model
1) measurement obtains
(1) asterism trajectory extraction and screening：The m frames that star sensor is obtained, m ＞ 2, blurred picture is handled, specifically Including removing ambient noise, all pixels point on every section of asterism track is extracted by connecting domain method, to all asterisms obtained Track is screened, and the situation that two or more asterism tracks overlap or overlap is avoided the occurrence of, simultaneously because the rotation of spacecraft Some asterisms may be caused to remove the visual field in imaging process, so needing to reject the point of image border, only retain whole expose Asterism track between light time all in visual field.
(2) asterism track fitting：By step, 1. obtained asterism track data is fitted and is obtained using Matlab programmings Fitting function is obtained, the direction vector of any time can be obtained by carrying out derivation as the fitting function to obtained by, and fuzzy to whole section Track carries out the size that integration obtains total Blur track, further calculates Optic flow information.
Assuming that kth frame star chart passes through the complete asterism track number that retains after image procossing, the i.e. asterism of complete exposure Number is n, then is considered as Z according to the Optic flow information of i-th star obtained after trajectory calculationi(k), it is expressed as：
Therefore, measurement Z (k) is represented by：
Z (k)=[Z1(k) Z2(k) ... Zn(k)] (9)
2) measurement model
System measurements model to establish process as follows：
Asterism light stream and the relation of three axis stabilized spacecraft angular speed are as follows：
Wherein, the coordinate that asterism is incident upon in image plane is (x, y), and f is the focal length of star sensor, and u and v are respectively x-axis , can be by being obtained to the derivation of asterism position with the light stream vector in y-axis direction.ω1, ω2And ω3Spacecraft is represented respectively to take the photograph Around x-axis under camera coordinate system, the angular speed of y-axis and z-axis.
According to the relation shown in formula (10), it is as follows to establish measurement model：
Z (k)=H (k) X (k)+V (k) (11)
H (k)=[H1(k)H2(k)…Hn(k)]T (12)
In formula, measure noise V (k) and be assumed to the white Gaussian noise that noise average is 0, xi(k) and yi(k) kth frame star is represented The imager coordinate of i-th star in figure.
3. Kalman filtering
Quantity of state is estimated by Kalman filter equation, filter is determined by the time for exposure of every frame star chart first Wave period T, is configured each parameter of Fast track surgery.Finally the initial parameter of Kalman filtering is determined.
By the use of formula (2) as state model, formula (11) is used as measurement model, with reference to the Optic flow information u (k) of extraction, v (k), spacecraft angular speed is estimated by Kalman filtering, finally obtains estimate and output information.
The principle of the present invention is：Using the information for obscuring star chart as foundation, multiframe is obscured after star chart carries out image procossing and carried The Optic flow information taken is established according to the relation between asterism light stream and three axis stabilized spacecraft angular speed as measurement and measures mould Type, state model selection Singer acceleration probabilistic models, three of spacecraft can be estimated finally by Kalman filtering The angular speed of axis.
The present invention compared with prior art the advantages of be：
(1) be directly based upon fuzzy star chart without carry out deblurring process, not against barycenter extraction and importance in star map recognition and Spacecraft angular speed is estimated with process, compensate for the deficiency of star sensor Attitude estimation under high-speed case, star can be made Sensor becomes the unique sensing element of spacecraft and provides measurement information without relying on gyro, reduces cost, avoids gyroscopic drift Influence, and its failure conditions can be prevented.
(2) with it is existing can only low speed smooth motion situation calculate optical flow approach compared with, can be in spacecraft high-speed motion Optic flow information is obtained in the fuzzy star chart of situation.
Brief description of the drawings
Fig. 1 is a kind of spacecraft Attitude rate estimator method flow diagram based on fuzzy star chart Optic flow information in the present invention.
Embodiment
Fig. 1 gives the flow chart of the spacecraft Attitude rate estimator method based on fuzzy star chart Optic flow information.In detail below Illustrate the specific implementation process of the present invention：
Fuzzy star chart of the spacecraft in the case of three different axis angular rates is obtained by star sensor first, this part with Constant rotational speed ω (t)=[0 °/s 0 °/s 5 °/s] of the spacecraft around three axisTExemplified by, since Kalman filtering is that a kind of recursion is calculated Method, to ensure its time renewal process, measures renewal process and its estimated accuracy, and fuzzy star chart frame number m must be more than 2 frames, so Selected simulation time is 30s, and the time for exposure per frame star chart is 1s, and 30 frames of generation obscure star chart.
1. the foundation of System State Model；
First, angular speed and angular acceleration of the quantity of state for spacecraft are selected
In formula, ω1(k), ω2(k) and ω3(k) when being kth frame star chart respectively spacecraft under camera coordinate system around x The angular speed of axis, y-axis and z-axis,WithRespectively kth frame star chart when spacecraft accelerate at the angle of three axis Degree.State model selects Singer acceleration probabilistic models, and state equation is as follows：
X (k)=Φ (T) X (k-1)+W (k-1) (2)
In formula, T is filtering cycle, Λ=diag { 1/ τ1,1/τ2,1/τ3, wherein τi(i=1,2,3), represents spacecraft The G-time constant of the corresponding x-axis of body, y-axis and z-axis, design parameter value are set in Kalman filtering process, and I is unit Matrix；W (k-1) is estimated system noise sequence suffered by state X (k), the variance matrix representation of W (k-1) for the kth frame star chart moment For Q (k-1), its each element is expressed as：
Q11(k-1)=Λ-2Σ2(4e-ΛT-3I-e-2ΛT+2ΛT) (4)
Q21(k-1)=Q12(k-1)=Λ-1Σ2(e-2ΛT+I-2e-ΛT) (5)
Q22(k-1)=Σ2(I-e-2ΛT) in (6) formula, Σ=diag { σ1,σ2,σ3} It is used for determining the population variance of noise in Singer acceleration probabilistic models, is expressed as：
Wherein, WithAngular acceleration element can be regarded asPdf model adjusting parameter, It is expressed as：Angular accelerationEqual to maximumAnd minimum valueProbability beProbability equal to 0 ForIn sectionOn be uniformly distributed, probability isParameter WithScope It is as follows：Occurrence is selected in Kalman filtering process.
2. measurement obtains and the foundation of system measurements model
1) measurement obtains
(1) asterism trajectory extraction and screening：The 30 frame blurred pictures that star sensor obtains are handled, are specifically included Except ambient noise, all pixels point on every section of asterism track is extracted by connecting domain method, to all asterism tracks for being obtained into Row screening, avoids the occurrence of the situation that two or more asterism tracks overlap or overlap, simultaneously because the rotation of spacecraft may be led Cause some asterisms to remove the visual field in imaging process, so needing to reject the point of image border, only retain the whole time for exposure Asterism track all in visual field.
(2) asterism track fitting：By step, 1. obtained asterism track data is fitted and is obtained using Matlab programmings Fitting function is obtained, the direction vector of any time can be obtained by carrying out derivation as the fitting function to obtained by, and fuzzy to whole section Track carries out the size that integration obtains total Blur track, further calculates Optic flow information.
Assuming that kth frame star chart passes through the complete asterism track number that retains after image procossing, the i.e. asterism of complete exposure Number is n, then is considered as Z according to the Optic flow information of i-th star obtained after trajectory calculationi(k), it is expressed as：
Therefore, measurement Z (k) is represented by：
Z (k)=[Z1(k)Z2(k)...Zn(k)] (9)
2) measurement model
System measurements model to establish process as follows：
Asterism light stream and the relation of three axis stabilized spacecraft angular speed are as follows：
Wherein, the coordinate that asterism is incident upon in image plane is (x, y), and f is the focal length of star sensor, and u and v are respectively x-axis , can be by being obtained to the derivation of asterism position with the light stream vector in y-axis direction.ω1, ω2And ω3Spacecraft is represented respectively to take the photograph Around x-axis under camera coordinate system, the angular speed of y-axis and z-axis.
According to the relation shown in formula (10), it is as follows to establish measurement model：
Z (k)=H (k) X (k)+V (k) (11)
H (k)=[H1(k)H2(k)…Hn(k)]T (12)
In formula, measure noise V (k) and be assumed to the white Gaussian noise that noise average is 0, xi(k) and yi(k) kth frame star is represented The imager coordinate of i-th star in figure.
3. Kalman filtering
Quantity of state is estimated by Kalman filter equation, filter is determined by the time for exposure of every frame star chart first Wave period T=1s, is arranged to each parameter of acceleration probabilistic model：
τi=10s,The initial parameter of Kalman filtering is defined below：
Original state X0For
X0=[0rad/s, 0rad/s, 0rad/s, 0rad/s2,0rad/s2,0rad/s2]T
Initial error variance matrix P0For：
P0=diag [10-5(rad/s)2,10-5(rad/s)2,10-5(rad/s)2,10-7(rad/s2)2,10-7(rad/s2 )2,10-7(rad/s2)2]
By the use of formula (2) as state model, formula (11) is used as measurement model, with reference to the Optic flow information u (k) of extraction, v (k), spacecraft angular speed is estimated by Kalman filtering, finally obtains estimate and output information.
The content not being described in detail in description of the invention belongs to the prior art known to professional and technical personnel in the field.

Claims (1)

  1. A kind of 1. spacecraft Attitude rate estimator method based on fuzzy star chart extraction Optic flow information, it is characterised in that：Select first It is the Singer acceleration probabilistic model of quantity of state as state model using spacecraft angular speed and angular acceleration, then to m frames Fuzzy star chart carries out image procossing, and m ＞ 2, the Optic flow information of extraction is as measurement, according to asterism light stream and spacecraft angular speed Between relation establish measurement model, the angular speed of three axis of spacecraft can be estimated finally by Kalman filtering, specific bag Include following steps：
    1) foundation of System State Model；
    First, angular speed and angular acceleration of the quantity of state for spacecraft are selected, i.e.,：
    <mrow> <mi>X</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>=</mo> <msup> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mrow> <msub> <mi>&amp;omega;</mi> <mn>1</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;omega;</mi> <mn>2</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <msub> <mi>&amp;omega;</mi> <mn>3</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <msub> <mover> <mi>&amp;omega;</mi> <mo>&amp;CenterDot;</mo> </mover> <mn>1</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <msub> <mover> <mi>&amp;omega;</mi> <mo>&amp;CenterDot;</mo> </mover> <mn>2</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mrow> <msub> <mover> <mi>&amp;omega;</mi> <mo>&amp;CenterDot;</mo> </mover> <mn>3</mn> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> </mfenced> <mi>T</mi> </msup> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>
    In formula, ω1(k), ω2(k) and ω3(k) it is kth frame respectively, 1≤k≤m, spacecraft is under camera coordinate system during star chart Around x-axis, the angular speed of y-axis and z-axis,WithRespectively kth frame star chart when spacecraft at the angle of three axis Acceleration；
    State model selects Singer acceleration probabilistic models, and state equation is as follows：
    X (k)=Φ (T) X (k-1)+W (k-1) (2)
    <mrow> <mi>&amp;Phi;</mi> <mrow> <mo>(</mo> <mi>T</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mi>I</mi> </mtd> <mtd> <mrow> <msup> <mi>&amp;Lambda;</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msup> <mrow> <mo>(</mo> <mi>I</mi> <mo>-</mo> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mi>&amp;Lambda;</mi> <mi>T</mi> </mrow> </msup> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mn>0</mn> </mtd> <mtd> <msup> <mi>e</mi> <mrow> <mo>-</mo> <mi>&amp;Lambda;</mi> <mi>T</mi> </mrow> </msup> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>3</mn> <mo>)</mo> </mrow> </mrow>
    In formula, T is filtering cycle, Λ=diag { 1/ τ1,1/τ2,1/τ3, wherein τi, i=1,2,3, represent spacecraft body phase The G-time constant of x-axis, y-axis and z-axis is answered, design parameter value is set in Kalman filtering process, and I is unit matrix；W (k-1) system noise sequence suffered by state X (k) is estimated for the kth frame star chart moment, the variance matrix of W (k-1) is expressed as Q (k- 1), its each element is expressed as：
    Q11(k-1)=Λ-2Σ2(4e-ΛT-3I-e-2ΛT+2ΛT) (4)
    Q21(k-1)=Q12(k-1)=Λ-1Σ2(e-2ΛT+I-2e-ΛT) (5)
    Q22(k-1)=Σ2(I-e-2ΛT) (6)
    In formula, Σ=diag { σ1,σ2,σ3It is the population variance for being used for determining noise in Singer acceleration probabilistic models, represent For：
    <mrow> <msubsup> <mi>&amp;sigma;</mi> <mi>i</mi> <mn>2</mn> </msubsup> <mo>=</mo> <mfrac> <msubsup> <mover> <mi>&amp;omega;</mi> <mo>&amp;CenterDot;</mo> </mover> <msub> <mi>M</mi> <mi>i</mi> </msub> <mn>2</mn> </msubsup> <mn>3</mn> </mfrac> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mn>4</mn> <msub> <mi>p</mi> <msub> <mi>M</mi> <mi>i</mi> </msub> </msub> <mo>-</mo> <msub> <mi>p</mi> <msub> <mn>0</mn> <mi>i</mi> </msub> </msub> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
    Wherein, WithAngular acceleration element can be regarded asPdf model adjusting parameter, i=1,2,3, expression For：Angular accelerationEqual to maximumAnd minimum valueProbability beI=1,2,3 is equal to 0 probabilityIn sectionOn be uniformly distributed, probability isParameter WithScope such as Under：Occurrence is selected in Kalman filtering process；
    2) measurement acquisition and the foundation of system measurements model
    (1) measurement obtains
    1. asterism trajectory extraction and screening：The m frames that star sensor is obtained, m ＞ 2, blurred picture is handled, and is specifically included Except ambient noise, all pixels point on every section of asterism track is extracted by connecting domain method, to all asterism tracks for being obtained into Row screening, avoids the occurrence of asterism track and overlaps or overlapping situation, due to the rotation of spacecraft may cause some asterisms into As process can remove the visual field, so needing to reject the point of image border, only retain star of the whole time for exposure all in visual field The locus of points；
    2. asterism track fitting：By step, 1. obtained asterism track data is fitted acquisition fitting function, by institute The fitting function obtained, which carries out derivation, can obtain the direction vector of any time, and integration is carried out to whole section of blurring trajectorie and obtains total mould The size in path is pasted, further calculates Optic flow information；
    Assuming that kth frame star chart is by the complete asterism track number that retains after image procossing, i.e. the asterism number of complete exposure N, then Z is considered as according to the Optic flow information of i-th star obtained after trajectory calculationi(k), it is expressed as：
    <mrow> <msub> <mi>Z</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mrow> <mi>u</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>v</mi> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>8</mn> <mo>)</mo> </mrow> </mrow>
    Therefore, measurement Z (k) is represented by：
    Z (k)=[Z1(k) Z2(k) ... Zn(k)] (9)
    (2) measurement model
    System measurements model to establish process as follows：
    Asterism light stream and the relation of three axis stabilized spacecraft angular speed are as follows：
    <mrow> <mtable> <mtr> <mtd> <mrow> <mi>u</mi> <mo>=</mo> <mi>f</mi> <mo>&amp;lsqb;</mo> <mo>-</mo> <mi>&amp;omega;</mi> <mfrac> <mrow> <mi>x</mi> <mi>y</mi> </mrow> <msup> <mi>f</mi> <mn>2</mn> </msup> </mfrac> <mo>+</mo> <msub> <mi>&amp;omega;</mi> <mn>2</mn> </msub> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mfrac> <msup> <mi>x</mi> <mn>2</mn> </msup> <msup> <mi>f</mi> <mn>2</mn> </msup> </mfrac> <mo>)</mo> </mrow> <mo>+</mo> <msub> <mi>&amp;omega;</mi> <mn>3</mn> </msub> <mfrac> <mi>y</mi> <mi>f</mi> </mfrac> <mo>&amp;rsqb;</mo> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mi>v</mi> <mo>=</mo> <mi>f</mi> <mo>&amp;lsqb;</mo> <mo>-</mo> <msub> <mi>&amp;omega;</mi> <mn>1</mn> </msub> <mrow> <mo>(</mo> <mn>1</mn> <mo>+</mo> <mfrac> <msup> <mi>y</mi> <mn>2</mn> </msup> <msup> <mi>f</mi> <mn>2</mn> </msup> </mfrac> <mo>)</mo> </mrow> <mo>+</mo> <msub> <mi>&amp;omega;</mi> <mn>2</mn> </msub> <mfrac> <mrow> <mi>x</mi> <mi>y</mi> </mrow> <msup> <mi>f</mi> <mn>2</mn> </msup> </mfrac> <mo>-</mo> <msub> <mi>&amp;omega;</mi> <mn>3</mn> </msub> <mfrac> <mi>x</mi> <mi>f</mi> </mfrac> <mo>&amp;rsqb;</mo> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>10</mn> <mo>)</mo> </mrow> </mrow>
    Wherein, the coordinate that asterism is incident upon in image plane is (x, y), and f is the focal length of star sensor, and u and v are respectively x-axis and y-axis The light stream vector in direction, can be by obtaining the derivation of asterism position；ω1, ω2And ω3Spacecraft is represented respectively to sit in video camera Around x-axis under mark system, the angular speed of y-axis and z-axis；
    According to the relation shown in formula (10), it is as follows to establish measurement model：
    Z (k)=H (k) X (k)+V (k) (11)
    H (k)=[H1(k) H2(k) … Hn(k)]T (12)
    <mrow> <msub> <mi>H</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfenced open = "[" close = "]"> <mtable> <mtr> <mtd> <mrow> <mo>-</mo> <mfrac> <mrow> <msub> <mi>x</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <msub> <mi>x</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> <mi>f</mi> </mfrac> </mrow> </mtd> <mtd> <mrow> <mi>f</mi> <mo>+</mo> <mfrac> <mrow> <msup> <msub> <mi>x</mi> <mi>i</mi> </msub> <mn>2</mn> </msup> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> <mi>f</mi> </mfrac> </mrow> </mtd> <mtd> <mrow> <msub> <mi>y</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> <mtr> <mtd> <mrow> <mo>-</mo> <mi>f</mi> <mo>-</mo> <mfrac> <mrow> <msup> <msub> <mi>y</mi> <mi>i</mi> </msub> <mn>2</mn> </msup> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> <mi>f</mi> </mfrac> </mrow> </mtd> <mtd> <mfrac> <mrow> <msub> <mi>x</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> <msub> <mi>y</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> <mi>f</mi> </mfrac> </mtd> <mtd> <mrow> <mo>-</mo> <msub> <mi>x</mi> <mi>i</mi> </msub> <mrow> <mo>(</mo> <mi>k</mi> <mo>)</mo> </mrow> </mrow> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> <mtd> <mn>0</mn> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>13</mn> <mo>)</mo> </mrow> </mrow>
    In formula, measure noise V (k) and be assumed to the white Gaussian noise that noise average is 0, xi(k) and yi(k) represent in kth frame star chart The imager coordinate of i-th star；
    3) Kalman filtering
    Quantity of state is estimated by Kalman filter equation, filtering week is determined by the time for exposure of every frame star chart first Phase T, is configured each parameter of Fast track surgery, finally the initial parameter of Kalman filtering is determined；
    Measurement model is used as by the use of formula (2) as state model, formula (9), with reference to the Optic flow information u (k) of extraction, v (k), Spacecraft angular speed is estimated by Kalman filtering, finally obtains estimate and output information.
CN201711456224.3A 2017-12-28 2017-12-28 A kind of spacecraft Attitude rate estimator method for extracting Optic flow information based on fuzzy star chart Active CN107942090B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201711456224.3A CN107942090B (en) 2017-12-28 2017-12-28 A kind of spacecraft Attitude rate estimator method for extracting Optic flow information based on fuzzy star chart

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201711456224.3A CN107942090B (en) 2017-12-28 2017-12-28 A kind of spacecraft Attitude rate estimator method for extracting Optic flow information based on fuzzy star chart

Publications (2)

Publication Number Publication Date
CN107942090A true CN107942090A (en) 2018-04-20
CN107942090B CN107942090B (en) 2019-10-29

Family

ID=61939523

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201711456224.3A Active CN107942090B (en) 2017-12-28 2017-12-28 A kind of spacecraft Attitude rate estimator method for extracting Optic flow information based on fuzzy star chart

Country Status (1)

Country Link
CN (1) CN107942090B (en)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108830807A (en) * 2018-06-01 2018-11-16 哈尔滨工业大学 A kind of star sensor image solution motion blur method based on MEMS gyro auxiliary
CN108955627A (en) * 2018-06-05 2018-12-07 北京航空航天大学 A method of the location information using sun sensor measurement carrier with respect to the sun
CN109029425A (en) * 2018-06-25 2018-12-18 中国科学院长春光学精密机械与物理研究所 A kind of fuzzy star chart restored method filtered using region
CN109343550A (en) * 2018-10-15 2019-02-15 北京航空航天大学 A kind of estimation method of the spacecraft angular speed based on moving horizon estimation

Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101082494A (en) * 2007-06-19 2007-12-05 北京航空航天大学 Self boundary marking method based on forecast filtering and UPF spacecraft shading device
CN101701822A (en) * 2009-11-06 2010-05-05 中国人民解放军国防科学技术大学 Star tracking method of star sensor based on correlation of optical joint and transformation
CN101846510A (en) * 2010-05-28 2010-09-29 北京航空航天大学 High-precision satellite attitude determination method based on star sensor and gyroscope
CN102435763A (en) * 2011-09-16 2012-05-02 中国人民解放军国防科学技术大学 Measuring method for attitude angular velocity of spacecraft based on star sensor
CN103148851A (en) * 2013-02-18 2013-06-12 清华大学 Method for determining attitude of star sensor based on roller shutter exposure imaging
CN103217544A (en) * 2013-03-21 2013-07-24 上海新跃仪表厂 Method and system for estimating star angular speed according to star point position change of star sensor
CN104280049A (en) * 2014-10-20 2015-01-14 北京控制工程研究所 Outfield precision testing method for high-precision star sensor
EP2180292B1 (en) * 2008-10-24 2015-02-18 JVC KENWOOD Corporation Apparatus and method for correcting the output signal of an angular velocity sensor
CN105374035A (en) * 2015-11-03 2016-03-02 北京航空航天大学 Star sensor star point extraction method under stray light interference
JP2017142344A (en) * 2016-02-09 2017-08-17 キヤノン株式会社 Image blur correction device, control method thereof, program, and storage medium

Patent Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101082494A (en) * 2007-06-19 2007-12-05 北京航空航天大学 Self boundary marking method based on forecast filtering and UPF spacecraft shading device
EP2180292B1 (en) * 2008-10-24 2015-02-18 JVC KENWOOD Corporation Apparatus and method for correcting the output signal of an angular velocity sensor
CN101701822A (en) * 2009-11-06 2010-05-05 中国人民解放军国防科学技术大学 Star tracking method of star sensor based on correlation of optical joint and transformation
CN101846510A (en) * 2010-05-28 2010-09-29 北京航空航天大学 High-precision satellite attitude determination method based on star sensor and gyroscope
CN102435763A (en) * 2011-09-16 2012-05-02 中国人民解放军国防科学技术大学 Measuring method for attitude angular velocity of spacecraft based on star sensor
CN103148851A (en) * 2013-02-18 2013-06-12 清华大学 Method for determining attitude of star sensor based on roller shutter exposure imaging
CN103217544A (en) * 2013-03-21 2013-07-24 上海新跃仪表厂 Method and system for estimating star angular speed according to star point position change of star sensor
CN104280049A (en) * 2014-10-20 2015-01-14 北京控制工程研究所 Outfield precision testing method for high-precision star sensor
CN105374035A (en) * 2015-11-03 2016-03-02 北京航空航天大学 Star sensor star point extraction method under stray light interference
JP2017142344A (en) * 2016-02-09 2017-08-17 キヤノン株式会社 Image blur correction device, control method thereof, program, and storage medium

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
袁洪琳 等: "基于可分离模糊核的复合运动模糊星图建模与仿真", 《红外与激光工程》 *

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108830807A (en) * 2018-06-01 2018-11-16 哈尔滨工业大学 A kind of star sensor image solution motion blur method based on MEMS gyro auxiliary
CN108830807B (en) * 2018-06-01 2022-01-28 哈尔滨工业大学 MEMS gyroscope-assisted star sensor image motion blur solving method
CN108955627A (en) * 2018-06-05 2018-12-07 北京航空航天大学 A method of the location information using sun sensor measurement carrier with respect to the sun
CN108955627B (en) * 2018-06-05 2020-09-11 北京航空航天大学 Method for measuring position information of carrier relative to sun by adopting sun sensor
CN109029425A (en) * 2018-06-25 2018-12-18 中国科学院长春光学精密机械与物理研究所 A kind of fuzzy star chart restored method filtered using region
CN109029425B (en) * 2018-06-25 2020-07-31 中国科学院长春光学精密机械与物理研究所 Fuzzy star map restoration method adopting regional filtering
CN109343550A (en) * 2018-10-15 2019-02-15 北京航空航天大学 A kind of estimation method of the spacecraft angular speed based on moving horizon estimation

Also Published As

Publication number Publication date
CN107942090B (en) 2019-10-29

Similar Documents

Publication Publication Date Title
US11057567B2 (en) Anti-shake method and apparatus for panoramic video, and portable terminal
CN107942090B (en) A kind of spacecraft Attitude rate estimator method for extracting Optic flow information based on fuzzy star chart
Rambach et al. Learning to fuse: A deep learning approach to visual-inertial camera pose estimation
CN111504312B (en) Unmanned aerial vehicle pose estimation method based on visual inertial polarized light fusion
CN106525074B (en) A kind of compensation method, device, holder and the unmanned plane of holder drift
CN105931275A (en) Monocular and IMU fused stable motion tracking method and device based on mobile terminal
Jia et al. Probabilistic 3-d motion estimation for rolling shutter video rectification from visual and inertial measurements
US11223764B2 (en) Method for determining bias in an inertial measurement unit of an image acquisition device
CN109540126A (en) A kind of inertia visual combination air navigation aid based on optical flow method
CN105953796A (en) Stable motion tracking method and stable motion tracking device based on integration of simple camera and IMU (inertial measurement unit) of smart cellphone
CN106814753B (en) Target position correction method, device and system
US20220051031A1 (en) Moving object tracking method and apparatus
CN108364319A (en) Scale determines method, apparatus, storage medium and equipment
CN107478220A (en) Unmanned plane indoor navigation method, device, unmanned plane and storage medium
CN107767425A (en) A kind of mobile terminal AR methods based on monocular vio
CN106709222B (en) IMU drift compensation method based on monocular vision
CN111609868A (en) Visual inertial odometer method based on improved optical flow method
CN105141807A (en) Video signal image processing method and device
CN105096337B (en) A kind of image global motion compensation method based on gyroscope hardware platform
CN109040525B (en) Image processing method, image processing device, computer readable medium and electronic equipment
CN111595332B (en) Full-environment positioning method integrating inertial technology and visual modeling
CN112204946A (en) Data processing method, device, movable platform and computer readable storage medium
CN112179373A (en) Measuring method of visual odometer and visual odometer
CN117073720A (en) Method and equipment for quick visual inertia calibration and initialization under weak environment and weak action control
CN114440877B (en) Asynchronous multi-camera visual inertial odometer positioning method

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