A method for calculating the position of satellite formation based on ranging information under the condition of star sensor correction

Abstract

The present invention discloses a method for solving the position of a satellite formation based on ranging information under the condition of star sensor correction: setting calibration points on the main star body and characteristic points on the slave star; solving the coordinates of the characteristic points of the slave star in the coordinate system of the main star, the relative position information of the slave star mass center to the main star mass center, the generalized vector of the position and velocity of the slave star mass center relative to the main star mass center in the main star coordinate system and its first-order derivative, entering the Kalman filter algorithm, and outputting the final satellite formation relative position solution result. The satellite formation position solution method based on ranging information under the condition of star sensor correction of the present invention takes laser ranging and attitude measurement as input, builds a state estimation equation based on a direct solution algorithm and relative kinematics equations, and uses the Kalman filter algorithm to give the relative position solution result of the satellite formation, thereby achieving high-precision inter-satellite relative navigation results of the satellite formation.

Description

Technical Field

[0001] The invention belongs to the technical field of satellite formation, and in particular relates to a satellite formation position solution method based on ranging information under the condition of star sensor correction. Background Art

[0002] Satellite formation technology plays an important role in the future aerospace field, and the relative position determination of satellite formations is the key guarantee for maintaining the formation configuration. Currently, most of the methods used are visual navigation methods, but this method will be affected when the distance is long, resulting in a decrease in imaging accuracy. At the same time, the distance error of visual imaging in the direction of the optical axis is large, which limits the use of visual methods in the process of determining the relative position of satellite formations. Compared with visual measurement, microwave ranging, laser ranging and other measurement methods have the characteristics of high precision, but their disadvantage is that ranging can usually only measure information in one dimension, while the relative position determination of satellite formations requires information in three dimensions. Therefore, a satellite formation position solution method that corrects the weak observability problem is urgently needed. Summary of the invention

[0003] In order to solve the above problems existing in the prior art, the present invention provides a satellite formation position solution method based on ranging information under the condition of star sensor correction. It takes laser ranging and attitude measurement (usually star sensor) as input, builds state estimation equation based on direct solution algorithm and relative kinematic equation, and uses Kalman filter algorithm to give the relative position solution result of satellite formation, realizes high-precision inter-satellite relative navigation result of satellite formation, and solves the problem of high-precision relative position determination of satellite formation.

[0004] The technical solution adopted by the present invention is a method for calculating the position of a satellite formation based on ranging information under the condition of star sensor correction, which specifically includes the following steps:

[0005] Step 1: Taking the coordinate system of the master star as the reference system, at least three master star calibration points with known coordinates are set on the master star, at least one slave star feature point is set on the slave star, and the distance between the slave star feature point and each master star calibration point is measured;

[0006] Step 2, according to the known coordinates of each main star calibration point in the main star coordinate system, and the distance between the slave star feature point and each main star calibration point, solve the coordinates of the slave star feature point in the main star coordinate system;

[0007] Step 3, according to the solved coordinates of the slave star feature point in the master star coordinate system and the known installation vector of the slave star feature point in the slave star coordinate system, solve the relative position information of the slave star center of mass to the master star center of mass;

[0008] Step 4, obtaining the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system and its first-order derivative;

[0009] Step 5: Based on the relative position information of the slave star mass center to the master star mass center obtained in step 3, and the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system obtained in step 4 and its first-order derivative, enter the Kalman filter algorithm and output the final satellite formation relative position solution result.

[0010] Furthermore, in step 1, the master star is equipped with three detector receiving ends, which serve as the three master star calibration points P1, P2, and P3, and their coordinates are P1(x1, y1, z1), P2(x2, y2, z2), and P3(x3, y3, z3), respectively; the slave star is equipped with a detector transmitting end, which serves as a slave star feature point Q. The coordinates of the slave star feature point Q in the master star coordinate system are Q(x Q ,y Q ,z Q ) The detector measures three distance information P1Q, P2Q, and P3Q, which are d1, d2, and d3 respectively.

[0011] Furthermore, in step 2, the method for solving the coordinates of the slave star feature point in the master star coordinate system specifically includes the following steps:

[0012] Define the generalized variable X to be solved Q =[x Q y Q z Q ] T ;X Q represents the solution of the coordinates of the slave star feature point Q in the master star coordinate system;

[0013] The initial coordinate value of the optional slave star feature point Q in the master star coordinate system

[0014] X Q Substitute the coordinates of the feature point of the slave star in the coordinate system of the master star into the solution function F(X Q ), as shown below:

[0015]

[0016] Perform i+1 iterations and convert the solution of the i-th iteration into Substitute the following formula:

[0017]

[0018] In the formula, It represents the solution result of the i+1th iteration, where i represents the number of iterations; Indicates that Substitute F(X Q ) The result obtained; Ai Equivalence Indicates the value of F(X Q ) to find X Q The Jacobian matrix of the first-order derivative of is calculated as follows:

[0019]

[0020] After multiple iterations, the coordinates of the slave star feature point Q in the master star coordinate system are obtained. Q .

[0021] Furthermore, in step 3, the relative position information from the star mass center to the master star mass center is solved as shown in the following formula:

[0022]

[0023] Where, O1 represents the coordinate of the mass center of the master star in the master star coordinate system; O2 represents the coordinate of the mass center of the slave star in the master star coordinate system; Represents the relative position information from the star mass center to the main star mass center; Represents the relative position information from the star feature point to the main star's center of mass; represents the installation vector of the slave satellite feature point Q in the slave satellite coordinate system; C1 and C2 represent the attitude transformation matrices of the master satellite and slave satellite coordinate systems relative to the inertial coordinate system respectively; the matrix operator T represents the transpose of the matrix.

[0024] Furthermore, in step 4, the generalized vectors of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system are as follows:

[0025]

[0026] The first-order derivative of the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system is as shown in the following formula:

[0027]

[0028] In the formula, X cw A generalized vector representing the position and velocity of the slave star's mass center relative to the master star's mass center in the master star's coordinate system; Represents X cw The first derivative of; x, y, z, They represent the position and velocity components of the slave star's mass center relative to the master star in the master star's coordinate system in three axes; They respectively represent the three-axis components of the acceleration of the center of mass of the slave star in the master star coordinate system; ω represents the instantaneous orbital angular velocity of the master star; ε represents the instantaneous orbital angular acceleration of the master star; μ represents the gravitational constant of the earth; and r represents the distance of the master star from the center of the earth.

[0029] Furthermore, in step 5, the Kalman filter algorithm is entered to output the final satellite formation relative position solution result, including the following steps:

[0030] Step 51, assigning the error variance matrix at the initial time: setting the generalized vector initial value X of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system cw (0), set the initial value of the error variance matrix P(0) = I6, I6 represents the sixth-order unit matrix; set the process noise matrix Q EKF =5×10 -5 I6; Set the measurement noise matrix R = 5 × 10 2 I6;

[0031] Step 52, enter the loop, and for the kth moment, calculate the state transfer matrix Φ(k):

[0032]

[0033]

[0034] Where Δt represents the filter step size, A represents the state transfer matrix, ω represents the instantaneous orbital angular velocity of the primary star; ε represents the instantaneous orbital angular acceleration of the primary star; μ represents the earth's gravitational constant; r represents the distance between the primary star and the center of the earth;

[0035] Step 53, calculate the one-step prediction at the kth moment:

[0036] X cw (k,k-1)=Φ(k)X cw (k-1)

[0037] P(k,k-1)=Φ(k)P(k-1)Φ(k)+Q EKF

[0038] In the formula, X cw (k,k-1) represents the kth moment X cw The system one-step prediction value, X cw (k-1) represents the k-1th moment X cw The system calculated value of; P(k,k-1) represents the one-step prediction value of the error variance matrix at the kth moment; P(k-1) represents the error variance matrix value at the k-1th moment;

[0039] Step 54, calculate the gain matrix K:

[0040] K=P(k,k-1)(P(k,k-1)+R) -1

[0041] Step 55, calculate the system measurement signal Z(k) at the kth moment:

[0042]

[0043] Step 56, error variance matrix update:

[0044] P(k)=(IK)P(k,k-1)(IK) T +KRK T

[0045] Where P(k) represents the error variance matrix value at the kth moment; I represents the 6th-order unit matrix;

[0046] Step 57, status update:

[0047] X cw (k) = X cw (k,k-1)+K[Z(k)-X cw (k,k-1)]

[0048] In the formula, X cw (k) represents the kth moment X cw The system value of

[0049] Step 58, the loop is completed at the kth moment, and the process enters the k+1th moment, or the loop is exited;

[0050] The Kalman filter algorithm is completed and the final satellite formation relative position solution output is output.

[0051] The beneficial effects of the embodiments of the present invention are as follows: the specific embodiments of the present invention use laser ranging and attitude measurement (usually star sensors) as input, build a state estimation equation based on a direct solution algorithm and relative kinematics equations, and use the Kalman filter algorithm to give a relative position solution result of a satellite formation, thereby achieving a high-precision inter-satellite relative navigation result of the satellite formation, solving the problem of the amplification effect of geometric characteristics on direct measurement errors and weak observability in the direct solution method, and providing a theoretical basis for the on-orbit operation of satellite formations. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0053] Figure 1 It is a schematic diagram of the measurement principle of the calibration point solution method in the prior art.

[0054] Figure 2It is a diagram of the error amplification effect of the calibration point solution method in the prior art.

[0055] Figure 3 It is a schematic diagram of solving the coordinates of the star mass center in the main star body coordinate system in an embodiment of the present invention.

[0056] Figure 4 It is a relative position estimation error curve of the satellite formation position solution method based on ranging information under the star sensor correction condition of the embodiment of the present invention. DETAILED DESCRIPTION

[0057] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0058] Each embodiment in this specification is described in a related manner, and the same or similar parts between the embodiments can be referred to each other, and each embodiment focuses on the differences from other embodiments. In particular, for the system embodiment, since it is basically similar to the method embodiment, the description is relatively simple, and the relevant parts can be referred to the partial description of the method embodiment.

[0059] The main existing solution method based on distance measurement information is the method of solving through calibration points, such as Figure 1 shown.

[0060] The satellite on the left is the master satellite, and the satellite on the right is the slave satellite. The transmitter of the detector carried by the slave satellite is point Q1 in the figure. The master satellite carries at least three detector receiving ends as shown in the figure P1, P2, and P3. The detectors can measure the three distance information P1Q1, P2Q1, and P3Q1. Considering that the three points P1P2P3 are known in the coordinate system of the master satellite, the coordinates of the slave satellite calibration point Q1 in the coordinate system of the master satellite can be solved. Considering the coordinates of Q1 in the coordinate system of the slave satellite and the attitude information of the two satellites, the relative position of the center of mass of the slave satellite relative to the master satellite can be obtained.

[0061] Theoretically, from the star calibration point Q1 to the jth calibration point P of the main star j (x j ,y j ,z j The distance of )(j=1,2,3) can be written as This quantity can also be measured by the ranging equipment (detector) carried by the slave and master satellites, that is, the measurement quantity d j , and then there is

[0062] It can be seen that in the above equation (x j ,y j ,z j ) and d j are all known quantities. The former is obtained when the satellite is assembled, and the latter is given by ranging information. Q ,y Q ,z Q ) are the coordinates to be solved. It can be seen that the above equation has three unknowns, so it is necessary to establish equations (Formula (1)) for the three calibration points and solve them.

[0063]

[0064] In the formula, (x Q ,y Q ,z Q ) represents the coordinates of the point Q1 to be solved, d j Indicates Q1 to P j The ranging information usually includes noise.

[0065] However, under some specific geometric configurations and motion conditions, such as when multiple satellite formations are in the same plane, the normal positioning accuracy of the plane will drop sharply due to the weak observability of the geometric configuration. For example, in the double-star formation mode, the initial value of the primary star is set in the geosynchronous orbit, with a rotation period around the earth of T1 = (23*60+56)*60+4s; the initial position of the primary star in the geocentric inertial coordinate system is r_t = (R1, 0, 0), where R1 = (T1·GM / 4π 2 )1 / 3 (M is the mass of the earth, G is the gravitational constant), the speed of the master star v_t=(0, V1, 0), where V1=(2π / T)*R1; the slave star T2=100*60s; the coordinate system of the master star is small, the position r_c=(R2, 0, 0), where R2=100m, the speed of the slave star v_c=(0, V2, 0), where V2=2πR 2 / T 2 , ranging accuracy 5×10 -3 m, and the coordinates of point Q1 are obtained by solving the calibration point equation through formula (1). The solution accuracy is as follows: Figure 2 As shown, it can be seen that the solution error is of the order of m at this time, and the ranging error is multiplied under specific geometric conditions (up to several thousand times), which is also a defect of the poor geometric precision factor of the calibration point calculation method.

[0066] Based on this, the present invention proposes a satellite formation position solution method based on ranging information under the condition of star sensor correction, and uses attitude measurement information to correct the weak observability problem, which specifically includes the following steps:

[0067] Step 1: Taking the coordinate system of the master star as the reference system, at least three master star calibration points with known coordinates are set on the master star, at least one slave star feature point is set on the slave star, and the distance between the slave star feature point and each master star calibration point is measured;

[0068] Step 2, according to the known coordinates of each main star calibration point in the main star coordinate system, and the distance between the slave star feature point and each main star calibration point, solve the coordinates of the slave star feature point in the main star coordinate system;

[0069] Step 3, according to the solved coordinates of the slave star feature point in the master star coordinate system and the known installation vector of the slave star feature point in the slave star coordinate system, solve the relative position information of the slave star center of mass to the master star center of mass;

[0070] Step 4, obtaining the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system and its first-order derivative;

[0071] Step 5: Based on the relative position information of the slave star mass center to the master star mass center obtained in step 3, and the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system obtained in step 4 and its first-order derivative, enter the Kalman filter algorithm and output the final satellite formation relative position solution result.

[0072] A specific embodiment of the present invention, as Figure 3 As shown in the figure, the satellite on the left is the master satellite, and the satellite on the right is the slave satellite. The master satellite carries three detector receivers, which serve as the three master star calibration points P1, P2, and P3. Their coordinates are known, namely P1(x1, y1, z1), P2(x2, y2, z2), and P3(x3, y3, z3). The slave satellite carries a detector transmitter, which serves as a slave star feature point Q. The coordinates of the slave star feature point Q in the master star coordinate system need to be solved, which is Q(x Q ,y Q ,z Q ). The detector measures and obtains three distance information P1Q, P2Q, and P3Q, which are d1, d2, and d3 respectively.

[0073] According to the known coordinates of the three main star calibration points P1, P2, P3 in the main star coordinate system, and the distances d1, d2, d3 between the slave star feature point and each calibration point, the coordinates of the slave star feature point Q in the main star coordinate system are solved. The specific solution method is:

[0074] First, define the generalized variable X Q =[x Q y Q z Q ] T ;X Qrepresents the coordinates of the slave star feature point Q in the master star coordinate system (x Q ,y Q ,z Q ) solution;

[0075] The initial coordinate value of the optional slave star feature point Q in the master star coordinate system Normally the initial value The selection can be determined by experience or the initial orbit parameters of the satellite. Usually, the closer the initial value is to the true value, the faster the iteration converges. In the specific embodiment of the present invention,

[0076] Then X Q Substitute the coordinates of the feature point of the slave star in the coordinate system of the master star into the solution function F(X Q ), F(X Q ) is the vector description of formula (1), as shown in formula (2):

[0077]

[0078] Perform i+1 iterations and convert the solution of the i-th iteration into Substitute into formula (3):

[0079]

[0080] In the formula, It represents the solution result of the i+1th iteration, where i represents the number of iterations; Indicates that Substitute F(X Q ) The result obtained; A i Equivalence Indicates the value of F(X Q ) to find X Q The Jacobian matrix of the first-order derivative of is calculated as shown in formula (4):

[0081]

[0082] Generally speaking, the greater the number of iterations, the higher the accuracy. In the specific embodiment of the present invention, the number of simulation iterations is selected as 50, that is, when the solution is The iteration ends when The output result is the solution X of the coordinates of the slave star feature point Q in the master star coordinate system. Q .

[0083] The coordinates of the master star mass center in the master star coordinate system are O1, and the coordinates of the slave star mass center in the master star coordinate system are O2. According to the solution of the coordinates of the slave star feature point Q in the master star coordinate system, Q, and the installation vector of the known slave star feature point Q in the slave star coordinate system, solve the relative position information of the slave star mass center to the master star mass center The specific calculation is shown in the following formula:

[0084]

[0085] In the formula, represents the relative position information from the star feature point to the center of mass of the master star, which is the solution of the coordinates of the slave star mark point in the master star coordinate system obtained by solving equations 2 to 4. Q ; represents the installation vector of the slave star feature point Q in the slave star coordinate system (since the vector does not change with the movement of the slave star, it can be obtained by installation and measurement on the ground before the satellite is launched, and is known information); C1 and C2 respectively represent the attitude conversion matrices of the master star and slave star coordinate systems relative to the inertial coordinate system (the specific embodiment of the present invention adopts the J2000 coordinate system), and the variables are measured by the star sensor (generally, the star sensor measures the attitude conversion matrix of its installation body relative to the inertial system, and its specific output forms include attitude quaternion, Euler angle, direction cosine matrix, etc., and the output of the direction cosine matrix is ​​adopted in the specific embodiment of the present invention); the matrix operator T represents the transpose of the matrix.

[0086] The kinematic law of the center of mass of the slave star relative to the center of mass of the master star conforms to the relative kinematic model of satellite relative motion, that is, the classic CW equation (Clohessy-Wiltshire) in satellite rendezvous and docking and orbit interception problems. The generalized vector of the position and velocity of the center of mass of the slave star relative to the center of mass of the master star in the coordinate system of the master star is X cw , as shown in formula (7), its first-order derivative Directly written as the CW equation form of equation (8):

[0087]

[0088]

[0089]

[0090] In the formula, Indicates the rate of change of the relative position of the slave star relative to the master star; Represents X cw The first derivative of cw The generalized vector representing the position and velocity of the center of mass of the slave star relative to the center of mass of the master star in the coordinate system of the master star; x, y, z, They represent the position and velocity components of the slave star's mass center relative to the master star in the master star's coordinate system in three axes; They respectively represent the three-axis components of the acceleration of the center of mass of the slave star in the master star coordinate system in the master star coordinate system; ω represents the instantaneous orbital angular velocity of the master star; ε represents the instantaneous orbital angular acceleration of the master star; μ represents the gravitational constant of the earth; r represents the distance of the master star from the center of the earth; the above three variables (ω, ε, r) are usually given by the orbit determination system (absolute navigation system) or the ground station.

[0091] According to the analytical solution of formula (8), we can get The specific analytical solution of formula (8) is common knowledge in the art.

[0092] The main purpose of this step is to: directly solve the model (2) ~ (4) under certain geometric conditions, there is a weak observation problem that will sharply amplify the ranging error, such as Figure 2 As shown, the micrometer-level ranging accuracy leads to a solution error of the meter level. In order to effectively utilize the high accuracy of laser ranging, the present invention uses high-precision attitude measurement (generally given by attitude measurement devices such as star sensors, that is, the attitude conversion matrices C1 and C2 in formula (5), which can be regarded as known matrices) to correct it, construct a kinematic model of the variables to be solved, and build a Kalman filter algorithm to process the measurement information, thereby suppressing the weak observability problem and providing a high-precision solution result.

[0093] Enter the Kalman filter algorithm to give the final high-accuracy position estimate X cw , including the following steps:

[0094] (1), the error variance matrix is ​​assigned at the initial moment:

[0095] Set the initial value of the generalized vector X of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system cw (0), X cw (0) is usually given by the integration of the inertial measurement element. This initial value allows for a large error. Set the initial value of the error variance matrix P(0) to be a sixth-order unit matrix, i.e., P(0) = I6. Set the process noise matrix Q EKF , for Q EKF =5×10 -5 I6; Set the measurement noise matrix to R = 5 × 10 2 I6;

[0096] (2) Enter the loop and calculate the state transfer matrix Φ(k) for the kth moment.

[0097]

[0098]

[0099] Where Δt represents the filtering step size, and A represents the state transfer matrix, which is obtained according to formula (9).

[0100] (3) Calculate the one-step prediction at the kth moment:

[0101]

[0102] In the formula, X cw (k,k-1) represents the kth moment X cw The system one-step prediction value of (can also be regarded as the transition value from time k-1 to time k), X cw (k-1) represents the k-1th moment X cw The system calculated value of; P(k,k-1) represents the one-step prediction value of the error variance matrix at the kth moment; P(k-1) represents the error variance matrix value at the k-1th moment;

[0103] (4), calculate the gain matrix K:

[0104] K=P(k,k-1)(P(k,k-1)+R) -1 (11)

[0105] (5), calculate the system measurement signal Z(k) at the kth moment:

[0106]

[0107] (6), the error variance matrix is ​​updated:

[0108] P(k)=(IK)P(k,k-1)(IK) T +KRK T (13)

[0109] Where P(k) represents the error variance matrix value at the kth moment; I represents the 6th-order identity matrix.

[0110] (7), Status update:

[0111] X cw (k) = X cw (k,k-1)+K[Z(k)-X cw (k,k-1)] (14)

[0112] In the formula, X cw (k) represents the kth moment X cw The system value of

[0113] (8) At the kth moment, the loop is completed and the process enters the k+1th moment, or the loop is exited (generally, the loop exit condition is receiving an instruction to stop the loop or the algorithm encounters a fault and crashes);

[0114] The Kalman filter algorithm is completed and the final satellite formation relative position solution result X is outputcw (k).

[0115] The relative position estimation accuracy is as follows: Figure 4 As shown in the figure, it can be seen that the relative position estimation accuracy of the method proposed in the present invention reaches 4×10 -3 m level, its estimation accuracy is comparable to that of the direct solution method ( Figure 2 ) is greatly improved. The satellite formation position calculation method based on ranging information under the star sensor correction condition of the present invention can effectively utilize the high precision of ranging means such as microwaves and lasers to improve the relative navigation accuracy, and at the same time solve the defect that ranging can only measure one dimension, and give a three-axis (x, y, z directions in the main star body coordinate system) relative navigation result.

[0116] The present invention directly models and solves (F(X) by calibrating points-feature points Q )) and relative kinematics modeling and solving (Equations 5, 6, and 7), establish a connection between the two, and correct the former with weak observation problems ( Figure 2 ). The specific embodiment of the present invention uses the relative dynamic equation of the slave star mass center relative to the master star mass center to establish a coupling relationship with the angle information measured by the star sensor, and further uses high-precision attitude measurement to correct the position solution result. If the calibration point-feature point directly solves the function F(X Q ) directly obtains the first-order derivative, which has a relatively complex form and structure, and it is difficult to introduce other high-precision measurement information as correction. Compared with the existing solution method, the solution of the present invention overcomes the problem of weak observability caused by special geometric configurations and greatly reduces the influence of geometric configurations on solution errors. This innovation can effectively ensure the high-precision formation maintenance of satellites and provide support for the development of future space missions (such as SAR formation imaging, optical interferometry imaging, gravitational wave detection, etc.).

[0117] The above description is only a preferred embodiment of the present invention and is not intended to limit the protection scope of the present invention. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention are included in the protection scope of the present invention.

Claims

1. A satellite formation position calculation method based on ranging information under star sensor correction conditions, characterized in that: The specific steps include: Step 1: Taking the coordinate system of the master star as the reference system, at least three master star calibration points with known coordinates are set on the master star, at least one slave star feature point is set on the slave star, and the distance between the slave star feature point and each master star calibration point is measured; Step 2, according to the known coordinates of each main star calibration point in the main star coordinate system, and the distance between the slave star feature point and each main star calibration point, solve the coordinates of the slave star feature point in the main star coordinate system; Step 3, according to the solved coordinates of the slave star feature point in the master star coordinate system and the known installation vector of the slave star feature point in the slave star coordinate system, solve the relative position information of the slave star center of mass to the master star center of mass; Step 4, obtaining the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system and its first-order derivative; Step 5: Enter the Kalman filter algorithm based on the relative position information of the slave star mass center to the master star mass center obtained in step 3, and the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system obtained in step 4, and its first-order derivative, and output the final satellite formation relative position solution result; In step 3, the relative position information from the centroid of the star to the centroid of the master star is solved as shown in the following formula: Where, O1 represents the coordinate of the mass center of the master star in the master star coordinate system; O2 represents the coordinate of the mass center of the slave star in the master star coordinate system; Represents the relative position information from the star mass center to the main star mass center; Represents the relative position information from the star feature point to the main star's center of mass; represents the installation vector of the slave satellite feature point Q in the slave satellite coordinate system; C1 and C2 represent the attitude transformation matrix of the master satellite and slave satellite coordinate systems relative to the inertial coordinate system respectively; the matrix operator T represents the transpose of the matrix; In step 4, the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system is as follows: The first-order derivative of the generalized vector of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system is as shown in the following formula: Where, X cw A generalized vector representing the position and velocity of the slave star's mass center relative to the master star's mass center in the master star's coordinate system; Represents X cw The first derivative of; x, y, z, They represent the position and velocity components of the slave star's mass center relative to the master star in the master star's coordinate system in three axes; They represent the three-axis components of the centroid acceleration of the slave star in the master star coordinate system; ω represents the instantaneous orbital angular velocity of the master star; ε represents the instantaneous orbital angular acceleration of the primary star; μ represents the gravitational constant of the Earth; r represents the distance between the primary star and the center of the Earth; In step 5, the Kalman filter algorithm is entered to output the final satellite formation relative position solution result, including the following steps: Step 51, assigning the error variance matrix at the initial time: setting the generalized vector initial value X of the position and velocity of the slave star mass center relative to the master star mass center in the master star coordinate system cw (0), set the initial value of the error variance matrix P(0) = I6, I6 represents the sixth-order unit matrix; set the process noise matrix Q EKF =5×10 -5 I6; Set the measurement noise matrix R = 5 × 10 2 I6; Step 52, enter the loop, and for the kth moment, calculate the state transfer matrix Φ(k): Where Δt represents the filter step size, A represents the state transfer matrix, ω represents the instantaneous orbital angular velocity of the primary star; ε represents the instantaneous orbital angular acceleration of the primary star; μ represents the earth's gravitational constant; r represents the distance between the primary star and the center of the earth; Step 53, calculate the one-step prediction at the kth moment: X cw (k,k-1)=Φ(k)X cw (k-1) P(k,k-1)=Φ(k)P(k-1)Φ(k)+Q EKF Where, X cw (k,k-1) represents the kth moment X cw The system one-step prediction value, X cw (k-1) represents the k-1th moment X cw The system calculated value of; P(k,k-1) represents the one-step prediction value of the error variance matrix at the kth moment; P(k-1) represents the error variance matrix value at the k-1th moment; Step 54, calculate the gain matrix K: K=P(k,k-1)(P(k,k-1)+R) -1 Step 55, calculate the system measurement signal Z(k) at the kth moment: Step 56, error variance matrix update: P(k)=(IK)P(k,k-1)(IK) T +KRK T Where P(k) represents the error variance matrix value at the kth moment; I represents the 6th-order unit matrix; Step 57, status update: X cw (k)=X cw (k,k-1)+K[Z(k)-X cw (k,k-1)] Where, X cw (k) represents the kth moment X cw The system value of Step 58, the loop is completed at the kth moment, and the process enters the k+1th moment, or the loop is exited; The Kalman filter algorithm is completed and the final satellite formation relative position solution output is output.

2. The satellite formation position calculation method based on ranging information under the star sensor correction condition according to claim 1 is characterized in that: In step 1, the master star is equipped with three detector receivers as the three master star calibration points P1, P2, and P3, whose coordinates are P1(x1, y1, z1), P2(x2, y2, z2), and P3(x3, y3, z3), respectively; the slave star is equipped with a detector transmitter as a slave star feature point Q. The coordinates of the slave star feature point Q in the master star coordinate system are Q(x Q ,y Q ,z Q ) The detector measures three distance information P1Q, P2Q, and P3Q, which are d1, d2, and d3 respectively.

3. The satellite formation position calculation method based on ranging information under the star sensor correction condition according to claim 1 is characterized in that: In step 2, the method for solving the coordinates of the slave star feature point in the master star coordinate system specifically includes the following steps: Define the generalized variable X to be solved Q =[x Q y Q z Q ] T ;X Q represents the solution of the coordinates of the slave star feature point Q in the master star coordinate system; The initial coordinate value of the optional slave star feature point Q in the master star coordinate system X Q Substitute the coordinates of the feature point of the slave star in the coordinate system of the master star into the solution function F(X Q ), as shown below: Perform i+1 iterations and convert the solution of the i-th iteration into Substitute the following formula: In the formula, It represents the solution result of the i+1th iteration, where i represents the number of iterations; Indicates that Substitute F(X Q ) The result obtained; A i Equivalence Indicates the value of F(X Q ) to find X Q The Jacobian matrix of the first-order derivative of is calculated as follows: After multiple iterations, the coordinates of the slave star feature point Q in the master star coordinate system are obtained. Q .

Images