CN107168350B - Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation - Google Patents

Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation Download PDF

Info

Publication number
CN107168350B
CN107168350B CN201710375156.1A CN201710375156A CN107168350B CN 107168350 B CN107168350 B CN 107168350B CN 201710375156 A CN201710375156 A CN 201710375156A CN 107168350 B CN107168350 B CN 107168350B
Authority
CN
China
Prior art keywords
angular velocity
target
frequency
service spacecraft
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.)
Expired - Fee Related
Application number
CN201710375156.1A
Other languages
Chinese (zh)
Other versions
CN107168350A (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.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical 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 Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN201710375156.1A priority Critical patent/CN107168350B/en
Publication of CN107168350A publication Critical patent/CN107168350A/en
Application granted granted Critical
Publication of CN107168350B publication Critical patent/CN107168350B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • BPERFORMING OPERATIONS; TRANSPORTING
    • B64AIRCRAFT; AVIATION; COSMONAUTICS
    • B64GCOSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00Cosmonautic vehicles
    • B64G1/22Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/244Spacecraft control systems
    • BPERFORMING OPERATIONS; TRANSPORTING
    • B64AIRCRAFT; AVIATION; COSMONAUTICS
    • B64GCOSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00Cosmonautic vehicles
    • B64G1/22Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/64Systems for coupling or separating cosmonautic vehicles or parts thereof, e.g. docking arrangements
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/14Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G06F17/141Discrete Fourier transforms
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Remote Sensing (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Discrete Mathematics (AREA)
  • Automation & Control Theory (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Computing Systems (AREA)
  • Navigation (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)

Abstract

The invention discloses a method for calculating the optimal rotation angular velocity of a service spacecraft during fixed-axis rotation, which comprises the following steps: firstly, carrying out Fourier transform on a measured attitude quaternion of a target to obtain three characteristic angular frequencies of target motion; then, a filtering method is used for enabling the estimation value of the angular frequency to be more accurate; and finally, multiplying the characteristic angular frequency with the maximum influence by 2 to obtain the optimal angular speed for the rotation of the service spacecraft. The self-rotation angular velocity designed according to the invention enables the service spacecraft to rotate, and can effectively reduce the relative angular velocity and the relative linear velocity between the service spacecraft and the target, thereby reducing the burden of mechanical arm in capturing and enhancing the safety of space on-orbit service. In addition, although the service spacecraft rotates on a fixed shaft, the self-rotation angular velocity is not limited to the target of the fixed shaft rotation, and the relative angular velocity can be effectively reduced for the target rolling along three shafts.

Description

Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation
The technical field is as follows:
the invention relates to an on-orbit service technology in the field of spaceflight, in particular to a method for calculating an optimal rotation angular speed of a service spacecraft during fixed-axis rotation.
Background art:
at present, the space on-orbit service technology is receiving general attention of researchers at home and abroad. Targets for space on-orbit services include faulty satellites and space debris. The aim is to repair the target using a service spacecraft or to clear the target out of the current orbit. One of the prerequisites for performing a spatial on-track service is to dock the targets. The external moment acting on the target in space is very small, so that the initial angular momentum of the target is difficult to dissipate, and the target is often in a rolling state. Docking of a tumbling object is particularly difficult because when the service spacecraft and the object have high relative speeds, collisions can easily occur to damage the docking mechanism. For a low-speed rolling satellite, the relative angular velocity can be eliminated by using the motion of the space manipulator, but when the target angular velocity is too large, the motion of the manipulator cannot meet the requirement. In the past, there has been a method of reducing the relative angular velocity of a target rotating on a fixed axis by rotating a service spacecraft body, but no method of eliminating the relative velocity to the maximum extent of a target rotating on a three-axis rolling motion has been proposed.
The invention content is as follows:
the invention aims to reduce the relative angular velocity of a service spacecraft and any rolling target, and provides a method for calculating the optimal rotation angular velocity of the service spacecraft during fixed-axis rotation, which is used for calculating the optimal rotation angular velocity of the service spacecraft during rotation in a space capture task.
In order to achieve the purpose, the invention adopts the following technical scheme to realize the purpose:
a method for calculating the optimal rotation angular velocity of a service spacecraft during fixed-axis rotation comprises the following steps:
1) measuring the real-time postures of the rolling target by using a visual observation system or a laser radar ranging system, recording by using a quaternion representation method, and simultaneously recording the time corresponding to each posture;
2) after acquiring observation information of N groups of target attitude quaternions, performing discrete Fourier transform on each component of the attitude quaternion, wherein due to the characteristics of an attitude kinetic equation, a plurality of peaks appear on the graph of the discrete Fourier transform result of the quaternion in a frequency domain, wherein the value of N needs to meet the requirement
Figure GDA0002196310620000021
fmSampling frequency, f, representing observed values in the time domain1The frequency corresponding to the first peak on the graph for transforming the observation data in the time domain into the frequency domain;
3) the frequencies corresponding to the three peaks with the highest peak height are used as initial values, and more quaternion observation data are utilized through an extended Kalman filtering method or other filtering methods, so that the estimated values of the frequencies are accurate;
4) after the relative error of the estimation of the frequency is less than 5%, finding the frequency corresponding to the peak with the highest peak value in the frequency domain graph of the quaternion, and directly calculating the optimal rotation angular velocity of the service spacecraft through the frequency.
The invention is further improved in that the observation system in the step 1) obtains a quaternion observation value q (t)k) Wherein t iskFor the time corresponding to the k-th observation, the observed value contains an error, and the error is considered as white gaussian noise for processing.
The invention is further improved in that the discrete Fourier transform method used in the step 2) is as follows:
Figure GDA0002196310620000022
wherein q isi(tk) Represents q (t)k) Z of the ith componenti(n) represents qi(tk) The value of the nth term in the frequency domain after discrete Fourier transform, N is the total number of observed values, j is an imaginary unit, zi(n) corresponding to a frequency in the frequency domain of
Figure GDA0002196310620000023
Wherein f ismThe sampling frequency representing the observed value in the time domain, with f (n) as the abscissa, zi(n) is the ordinate, i.e. q can be plottediIn the frequency domain, the pattern shows several peaks, wherein the discrete Fourier transform is such that
Figure GDA0002196310620000031
The abscissa, i.e. frequency, corresponding to the first three highest peaks is respectively denoted as f1,f2And f3Wherein f is1The corresponding peak is highest.
A further development of the invention is that in the filtering method used in step 3), the quantity of states in the filtering process is x ═ uTT]TWherein u ═ u1,u2,…,u24]TIs aConstant vector, θ ═ θ123]TFor the vector composed of characteristic angular frequencies, before the filter is operated, the initial value of u is a zero vector, and the initial value of theta is theta0=2π[f1,f2,f3]T(ii) a The initial value of the error variance matrix of the state quantity is
Figure GDA0002196310620000032
Wherein I24×24And I3×3The unit matrixes with the dimension of 24 multiplied by 24 and the dimension of 3 multiplied by 3 respectively, and the one-step prediction equation of the state vector is
xk=xk-1
Wherein xkAnd xk-1Respectively, the state quantity x is at tkTime t andk-1the estimated value of the time, and since the quaternion is observable, the observation equation based on the state vector is
Figure GDA0002196310620000033
During the filtering process, theta1The estimated value of (c) is more and more accurate.
The invention further improves that the optimal rotation angular speed of the spacecraft finally obtained in the step 4) is 2 theta1And the direction of the angular velocity vector coincides with the angular momentum direction of the target, and the service spacecraft approaches to and catches the target along the angular momentum direction of the target.
The invention has the following beneficial effects:
according to the self-rotation angular velocity designed by the invention, the service spacecraft rotates, the relative angular velocity and the relative linear velocity between the service spacecraft and the target are reduced, the burden of a mechanical arm during capturing can be reduced, and the safety of space on-orbit service is enhanced. Because the service spacecraft rotates at constant speed with the fixed shaft, extra torque is not required to be consumed to control the rotating speed in the capturing process, and the energy consumption on the service spacecraft can be saved. In addition, the optimal rotating speed is calculated only by using a sensor of the existing service spacecraft, and extra load is not required to be added on the service spacecraft.
Furthermore, the invention allows the existence of observation errors in the calculation process of the optimal rotating speed, and is more suitable for the situation in practical application.
Furthermore, the invention firstly uses a Fourier transform method to determine a rough value of the optimal angular frequency, which can be an initial value of the extended Kalman filtering method, is favorable for convergence of the extended Kalman filtering method, and enhances the stability of the invention in practical application.
Furthermore, the invention uses the extended Kalman filter to calculate the optimal rotating speed in real time, thereby not only eliminating the influence of observation errors, but also ensuring that the calculation precision of the optimal rotating speed is more and more accurate along with the increase of the observed quantity, so that the rotating speed control of the service spacecraft is not needed to be carried out after the calculation is finished, and the optimal rotating speed can be calculated while being controlled until the optimal rotating speed is more and more close to the optimal rotating speed.
In addition, although the service spacecraft rotates on a fixed axis, because the invention accurately estimates the most obvious angular frequency component when the rolling target rotates, the invention designs the self-rotation angular velocity to be not only suitable for the target rotating on the fixed axis, but also effectively reduce the relative angular velocity for the target rolling along three axes.
Description of the drawings:
FIG. 1: schematic diagram of docking with a tumbling target using a rotating service spacecraft.
FIG. 2: and the corresponding relation between the Fourier transform result of the attitude quaternion of the rolling target and the angular frequency.
FIG. 3: and (b) comparing the projection of the target butt joint shaft end points in the body coordinate system of the service spacecraft when the service spacecraft is static (a) and rotates at the optimal rotating speed (b).
FIG. 4: and (b) comparing the end points of the butt joint shaft of the target with the motion speed of the service spacecraft when the service spacecraft is static (a) and rotates at the optimal rotating speed (b).
The specific implementation mode is as follows:
the invention is further illustrated by the following figures and examples.
FIG. 1 is a schematic illustration of docking with a tumbling target using a rotating service spacecraft. As shown, the serving spacecraft rotates on a fixed axis, so that its angular momentum coincides with the angular velocity direction and points to the angular momentum direction of the target. Since the target approximation is not considered to be affected by external moments, its angular momentum is conserved, and the orientation of the angular momentum direction in the inertial frame is fixed. Meanwhile, since the target makes an arbitrary rolling motion, its angular velocity vector does not necessarily coincide with the angular momentum direction, and its orientation in the inertial coordinate system is generally changed with time.
Through a non-contact measurement method, such as a binocular vision method, a laser ranging radar method and the like, the service spacecraft can measure the change condition of the attitude quaternion of the target along with the time. In this example, the initial state of the target is taken to be
Figure GDA0002196310620000051
ω0=[4 5 6]Trad/s
q0=[1 0 0 0]T
Wherein I is the representation of the inertia tensor of the target in the principal axis coordinate system of the body, omega0Representation of the target initial angular velocity in the principal axis coordinate system of the body, q0The initial value of the attitude quaternion of the body principal axis coordinate system relative to the inertial coordinate system is the target.
The observed quaternion data is fourier-transformed, and the frequency corresponding to the abscissa in the frequency domain is multiplied by 2 pi, so that an image of the fourier-transformed result corresponding to the angular frequency can be obtained, as shown in fig. 2. The angular frequencies corresponding to the peaks in the graph are the three characteristic angular frequencies.
Further obtaining the accurate value of the angular frequency with the highest peak value by using a filtering method to obtain theta13.6583 rad/s. Then obtaining the optimal rotating speed of the service spacecraft as 2 theta17.3166 rad/s. Assuming the target docking axis is in the z-direction of the ontology coordinate axis, FIG. 3 shows the service spacecraft at rest (a)And (b) projecting the target butt joint shaft end point in the body coordinate system of the service spacecraft when rotating at the optimal rotating speed. It can be seen that the motion trajectory of the end point of the docking shaft is simplified by autorotation. Meanwhile, fig. 4 shows the moving speed of the butt-joint shaft end point relative to the service spacecraft, which is targeted when the service spacecraft is stationary (a) and when it spins at the optimal rotation speed (b). It can be seen that by spinning, the relative speed of movement of the end points of the stub shafts is much reduced. In summary, the service spacecraft can greatly reduce the butt joint difficulty of the target spindle end point through autorotation, thereby improving the success rate of butt joint.

Claims (5)

1. A method for calculating the optimal rotation angular velocity of a service spacecraft during fixed-axis rotation is characterized by comprising the following steps:
1) measuring the real-time postures of the rolling target by using a visual observation system or a laser radar ranging system, recording by using a quaternion representation method, and simultaneously recording the time corresponding to each posture;
2) after acquiring observation information of N groups of target attitude quaternions, performing discrete Fourier transform on each component of the attitude quaternion, wherein due to the characteristics of an attitude kinetic equation, a plurality of peaks appear on the graph of the discrete Fourier transform result of the quaternion in a frequency domain, wherein the value of N needs to meet the requirement
Figure FDA0002372253040000011
fmSampling frequency, f, representing observed values in the time domain1The frequency corresponding to the first peak on the graph for transforming the observation data in the time domain into the frequency domain;
3) the frequencies corresponding to the three peaks with the highest peak height are used as initial values, and more quaternion observation data are utilized through an extended Kalman filtering method or other filtering methods, so that the estimated values of the frequencies are accurate;
4) after the relative error of the estimation of the frequency is less than 5%, finding the frequency corresponding to the peak with the highest peak value in the frequency domain graph of the quaternion, and directly calculating the optimal rotation angular velocity of the service spacecraft through the frequency.
2. The method for calculating the optimal rotation angular velocity of a serving spacecraft during fixed-axis rotation according to claim 1, wherein the observation system in step 1) obtains a quaternion observation value q (t)k) Wherein t iskFor the time corresponding to the k-th observation, the observed value contains an error, and the error is considered as white gaussian noise for processing.
3. The calculation method for the optimal rotation angular velocity of the serving spacecraft during the fixed axis rotation according to claim 2, wherein the discrete fourier transform method used in the step 2) is as follows:
Figure FDA0002372253040000012
wherein N is 1,2, …, N, i is 0,1,2,3, qi(tk) Represents q (t)k) Z of the ith componenti(n) represents qi(tk) The value of the nth term in the frequency domain after discrete Fourier transform, N is the total number of observed values, j is an imaginary unit, zi(n) corresponding to a frequency in the frequency domain of
Figure FDA0002372253040000021
Wherein f ismThe sampling frequency representing the observed value in the time domain, with f (n) as the abscissa, zi(n) is the ordinate, i.e. q can be plottediIn the frequency domain, the pattern shows several peaks, wherein the discrete Fourier transform is such that
Figure FDA0002372253040000022
The abscissa, i.e. frequency, corresponding to the first three highest peaks is respectively denoted as f1,f2And f3Wherein f is1The corresponding peak is highest.
4. A method as claimed in claim 3, wherein the filtering method used in step 3) is such that the state quantity during the filtering process is x ═ uTT]TWherein u ═ u1,u2,…,u24]TIs a vector of constants, [ theta ] - [ theta ]123]TFor the vector composed of characteristic angular frequencies, before the filter is operated, the initial value of u is a zero vector, and the initial value of theta is theta0=2π[f1,f2,f3]T(ii) a The initial value of the error variance matrix of the state quantity is
Figure FDA0002372253040000023
Wherein I24×24And I3×3The unit matrixes with the dimension of 24 multiplied by 24 and the dimension of 3 multiplied by 3 respectively, and the one-step prediction equation of the state vector is
xk=xk-1
Wherein xkAnd xk-1Respectively, the state quantity x is at tkTime t andk-1the estimated value of the time, and since the quaternion is observable, the observation equation based on the state vector is
Figure FDA0002372253040000024
During the filtering process, theta1The estimated value of (c) is more and more accurate.
5. The method for calculating the optimal rotation angular velocity of the service spacecraft in the process of fixed-axis rotation according to claim 4, wherein the optimal rotation angular velocity of the spacecraft finally obtained in the step 4) is 2 θ1And the direction of the angular velocity vector coincides with the angular momentum direction of the target, and the service spacecraft approaches to and catches the target along the angular momentum direction of the target.
CN201710375156.1A 2017-05-24 2017-05-24 Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation Expired - Fee Related CN107168350B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710375156.1A CN107168350B (en) 2017-05-24 2017-05-24 Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710375156.1A CN107168350B (en) 2017-05-24 2017-05-24 Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation

Publications (2)

Publication Number Publication Date
CN107168350A CN107168350A (en) 2017-09-15
CN107168350B true CN107168350B (en) 2020-04-21

Family

ID=59821836

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710375156.1A Expired - Fee Related CN107168350B (en) 2017-05-24 2017-05-24 Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation

Country Status (1)

Country Link
CN (1) CN107168350B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111288888B (en) * 2018-12-10 2021-08-10 中国科学院沈阳自动化研究所 Large-size circular ring target structured light measuring method for automatic capture by manipulator

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP5569681B2 (en) * 2010-04-23 2014-08-13 国立大学法人 東京大学 POSITION ESTIMATION DEVICE AND POSITION ESTIMATION METHOD OF MOBILE BODY USING INDUCTION SENSOR, MAGNETIC SENSOR AND SPEED METER
US9073648B2 (en) * 2013-02-15 2015-07-07 The Boeing Company Star tracker rate estimation with kalman filter enhancement
CN103218482B (en) * 2013-03-29 2017-07-07 南京航空航天大学 The method of estimation of uncertain parameter in a kind of dynamic system
CN103983278B (en) * 2014-05-19 2016-07-27 中国人民解放军国防科学技术大学 A kind of measure the method affecting Satellite Attitude Determination System precision
CN106482896A (en) * 2016-09-28 2017-03-08 西北工业大学 A kind of contactless factor of inertia discrimination method of arbitrary shape rolling satellite
CN106468554B (en) * 2016-09-29 2018-05-15 西北工业大学 A kind of measuring method of the inertial parameter of contactless rolling satellite

Also Published As

Publication number Publication date
CN107168350A (en) 2017-09-15

Similar Documents

Publication Publication Date Title
CN108227492B (en) Identification method for tail end load dynamic parameters of six-degree-of-freedom series robot
CN107607947B (en) On-line estimation method for imaging parameters of satellite-borne radar based on Kalman filtering
CN105976353A (en) Spatial non-cooperative target pose estimation method based on model and point cloud global matching
CN107421541B (en) Method for measuring and calculating morphological parameters of fault-tolerant non-contact failure satellite
CN110125936A (en) A kind of the Shared control method and ground experiment verifying system of robot for space
CN103217544B (en) Method and system for estimating star angular speed according to star point position change of star sensor
CN108680198B (en) Relative navigation target inertia parameter identification method based on plume disturbance
CN109426147B (en) Adaptive gain adjustment control method for combined spacecraft after satellite acquisition
CN110503713B (en) Rotation axis estimation method based on combination of trajectory plane normal vector and circle center
CN110308459B (en) Model-independent non-cooperative satellite relative pose measurement method
Abdelrahman et al. Sigma-point Kalman filtering for spacecraft attitude and rate estimation using magnetometer measurements
CN107702709A (en) A kind of noncooperative target moves the time-frequency domain mixing discrimination method with inertial parameter
CN107167145B (en) Form parameter measuring and calculating method of self-adaptive non-contact failure satellite
CN107246883A (en) A kind of Rotating Platform for High Precision Star Sensor installs the in-orbit real-time calibration method of matrix
CN109093620B (en) Binocular camera assisted space non-cooperative target kinetic parameter identification method
CN106482896A (en) A kind of contactless factor of inertia discrimination method of arbitrary shape rolling satellite
CN110567462B (en) Identification method for three-axis rotational inertia ratio of approximate spinning non-cooperative spacecraft
CN107168350B (en) Calculation method for optimal rotation angular velocity of service spacecraft during fixed-axis rotation
Liu et al. Mass and mass center identification of target satellite after rendezvous and docking
Biondi et al. Kinematic registration and shape analysis for locating center of mass in large passive spacecraft
Feng et al. Pose and motion estimation of unknown tumbling spacecraft using stereoscopic vision
EP2879011A1 (en) On-board estimation of the nadir attitude of an Earth orbiting spacecraft
CN111412919A (en) Method and device for calculating initial orbit error of space target
CN109145387B (en) Intelligent identification method of space rolling target inertia characteristics based on characteristic frequency
Xu et al. Vision-based moment of inertia estimation of non-cooperative space object

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
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20200421