A laplace shortest arc orbit determination algorithm based on particle swarm optimization

Abstract

The application discloses a Laplace extremely short arc orbit determination algorithm based on a particle swarm algorithm, and relates to the field of short arc orbit determination calculation of a spacecraft. The calculation method comprises the following steps: acquiring a coordinate position vector of an observed satellite at each moment and a direction vector obtained by observation; calculating a Kepler parameter initial value corresponding to the direction vector at each moment and a target position vector corresponding to the direction vector in an inertial system by using a particle swarm algorithm; performing orbit determination, correction and outlier elimination on the target position vector at each moment by using a Laplace method and a least square algorithm; and calculating orbit parameters corresponding to a short arc segment and a true anomaly corresponding to each moment, so that orbit determination calculation from the direction vector to the final orbit parameter of the target is realized, deviation caused by outliers during observation is offset, and extremely short arc orbit determination can be rapidly performed.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of spacecraft short-arc orbit determination, and particularly relates to a Laplace short-arc orbit determination algorithm based on a particle swarm algorithm. BACKGROUND

[0002] When a space-based detection satellite uses a camera to observe a space target, the space target is photographed with a star as a background. Since the direction of the star is known, the relative positions of the space target and the star on an image can be measured by a coordinate measuring instrument to obtain the azimuth information of the space target. Although an optical payload can be used to detect a space target to obtain an angle measurement, the distance between the detector and the debris cannot be obtained, and the target needs to be orbitally determined and cataloged based on observation data within a short observation time, which puts certain requirements on the input of the algorithm, that is, the target needs to be orbitally determined for a very short arc in a limited number of observation directions.

[0003] At present, there are different orbit determination methods for different numbers of stations (single, multiple) and different properties of measurement quantities (position, velocity, angle). For very short-arc orbit determination, common methods include Kepler element method, differential correction method, optimization algorithm and extended Kalman filter algorithm. The Kepler element method is simple in principle and fast in calculation, but has low accuracy and is difficult to remove outliers. The differential correction method and the optimization algorithm have the disadvantage of high computational complexity, and the extended Kalman filter algorithm requires a high accuracy of the original dynamic model.

[0004] Therefore, there is an urgent need for a Laplace short-arc orbit determination algorithm based on a particle swarm algorithm to solve the above problems. SUMMARY

[0005] To solve the above problems, the present application provides a Laplace short-arc orbit determination algorithm based on a particle swarm algorithm, which aims to quickly determine a very short arc and eliminate the influence of outliers and residuals to achieve accurate orbit determination.

[0006] In the present application, a Laplace short-arc orbit determination algorithm based on a particle swarm algorithm is provided, which includes the following steps:

[0007] S1, obtaining the coordinate position vector of the observation satellite at each time and the direction vector obtained by observation, and using a particle swarm algorithm to calculate the Kepler parameter initial value (a, e, M) corresponding to the direction vector at each time and the target position vector corresponding to the direction vector in the inertial system;

[0008] S2, using the Laplace method and the least squares algorithm to orbitally determine, correct and eliminate outliers of the target position vector obtained in S1 at each time, and calculating the orbit parameters corresponding to the short-arc segment and the true anomaly corresponding to each time.

[0009] Preferably, the specific steps of S1 are as follows:

[0010] S101: Obtain three variables in the range of Kepler parameters to obtain the particle swarm algorithm variables and define the three variables as the initial values of the Kepler parameters (a, e, M);

[0011] S102: Calculate the plane angle and the corresponding true angle at each time according to the Kepler parameter values, and calculate the position vector module length of the target according to the plane angle and the corresponding true angle;

[0012] S103: Obtain the direction vector obtained by observation and the coordinate position vector of the observation satellite, and calculate the included angle between the two to obtain the initial estimated value of the position vector of the target satellite;

[0013] S104: Calculate the target fitness function corresponding to the particle swarm algorithm according to the initial estimated value of the position vector obtained in S103, and update the particle swarm algorithm variables;

[0014] S105: Repeat the operations of S101 to S103 to perform iteration, obtain the optimal target fitness function value, and output the initial value of the Kepler parameter and the position vector value at each time.

[0015] Preferably, in step S102, the calculation expressions of the plane angle and the module length are as follows:

[0016] M i =M0+n(t i -t0);

[0017] |R i |=a(1-ecosf i );

[0018] Wherein, M is the plane angle value, n is the average angular velocity of the orbit, a is the semi-major axis value of the orbit, e is the eccentricity, R is the target position vector, is the true angle, i is the position angle index.

[0019] Preferably, the expression used in S103 to calculate the position vector value is as follows:

[0020]

[0021] R i =ρ i ×nj i +r i ;

[0022] Wherein, r i is the coordinate vector of the observation satellite in the earth-centered inertial system, nj i is the direction vector obtained by observation for the target, β is the included angle between the two vectors, ρ is the estimated value of the position vector in the Kepler parameter i ρ is the estimated value of the vector length from the observation satellite to the target satellite.

[0023] Preferably, the expression of the target fitness function in S104 is as follows:

[0024]

[0025] Wherein, R is the calculated position vector, F is the particle swarm fitness value corresponding to each iteration, and N is the number of observation times.

[0026] Preferably, the specific steps of S2 are as follows:

[0027] S201: Adopting polynomial curve smoothing fitting, the position vector values at several times are corrected;

[0028] S202: Obtaining the position and velocity information at t0, and calculating the initial values of the orbit semi-major axis, the argument of perigee, the eccentric anomaly, the parameter p and the parameter ;

[0029] S203: Using the values obtained in S202, the Laplace function value is calculated, and the orbit is extrapolated to obtain the position vector values and velocity vector values at several times;

[0030] S204: The position vector values at several times are substituted into S201, and the contents of S201 to S203 are repeated for iterative calculation to obtain the accurate value and the accurate value;

[0031] S205: Using the position vector and velocity vector values, the orbit elements are calculated and the orbit is determined.

[0032] Preferably, the specific content of adopting polynomial curve smoothing fitting in S201 is as follows:

[0033] Obtaining the initial position vector and the velocity vector at the time t0;

[0034] Adopting polynomial fitting function for fitting;

[0035] The expression of the polynomial fitting function is as follows:

[0036]

[0037] Wherein, R is the position vector value, t is the time value, d is the position vector at each time, j is the position index,

[0038] Preferably, the expression of the Laplace function adopted in S203 is as follows:

[0039]

[0040] Wherein, V is the speed of the target, f, g, F, G are Laplace operators.

[0041] Preferably, the orbit extrapolation expression in S203 is as follows:

[0042]

[0043] Wherein, a is the semi-major axis of the orbit, r0 is the position vector module length value, v0 is the velocity vector module length value, and p is the Laplace vector module length value.

[0044] Preferably, in the iterative calculation process in S204, when and the iteration is terminated.

[0045] In summary, compared with the traditional spacecraft short-arc orbit determination calculation technology, the Laplace extremely short-arc orbit determination algorithm based on the particle swarm algorithm of the present application is composed of the particle swarm algorithm, the Laplace method of orbit extrapolation and the least square algorithm, the preliminary operation of the direction vector to the Kepler constant and the position information at each time is completed by the particle swarm algorithm, the influence of the wild value and the residual is eliminated by the orbit data correction of the observation vector through the Laplace method and the least square method, and the relatively accurate orbit parameter value is obtained through the mutual iteration between the two, the traditional particle swarm algorithm orbit determination calculation is improved, the accuracy of the multi-vector orbit determination calculation is improved, and good applicability is obtained.

[0046] The technical solutions of the present application will be further described in detail below through the drawings and examples. BRIEF DESCRIPTION OF DRAWINGS

[0047] Figure 1 The flowchart of the Laplace extremely short-arc orbit determination algorithm based on the particle swarm algorithm of the present application is shown in the figure.

[0048] Figure 2 The orbit calculation example result graph of the Laplace extremely short-arc orbit determination algorithm based on the particle swarm algorithm of the present application is shown in the figure. DETAILED DESCRIPTION

[0049] The technical solutions of the present application will be further described in detail below through the drawings and examples.

[0050] The following description of at least one exemplary embodiment is merely illustrative in nature and is in no way limiting to the scope of the application or its applications or uses.

[0051] Techniques, methods, and apparatus known to those of ordinary skill in the relevant art can have not been discussed in any detail, but should be considered as part of the specification where appropriate.

[0052] In all of the compositions shown and discussed herein, any specific values should be interpreted as merely exemplary, and not a limitation. Thus, other examples of the exemplary embodiments can have different values.

[0053] Unless otherwise defined, technical terms or scientific terms used in the present application shall have the ordinary meanings as understood by one of ordinary skill in the art to which the present application pertains.

[0054] As shown in Figure 1 The present application provides a Laplace short arc orbit determination algorithm based on particle swarm algorithm, specifically as follows, S1, obtaining the coordinate position vector of each moment of observation satellite and the direction vector obtained by observation, using particle swarm algorithm to calculate the corresponding Kepler parameter initial value (a, e, M) of each moment of direction vector and the corresponding target position vector in the inertial system.

[0055] Further, the specific steps of S1 are as follows, S101: obtaining three variables in the range of Kepler parameters to obtain particle swarm algorithm variables and defining the three variables as the initial values of Kepler parameters (a, e, M).

[0056] S102: calculating the plane angle and the corresponding true angle of each moment according to the Kepler parameter value, and calculating the length of the position vector of the target according to the plane angle and the corresponding true angle.

[0057] Further, in step S102, the expressions for calculating the plane angle and the length are as follows:

[0058] M i = M0+ n(t i -t0).

[0059] |R i | = a (1-ecosf i ).

[0060] Wherein, M is the plane angle value, n is the average angular velocity of the orbit, a is the semi-major axis value of the orbit, e is the eccentricity, R is the target position vector, is the true angle, and i is the position angle index.

[0061] S103: obtaining the direction vector obtained by observation and the coordinate position vector of the observation satellite and calculating the included angle between the two to obtain the initial estimated value of the position vector of the target satellite.

[0062] Further, the expression used by S103 to calculate the position vector value is as follows:

[0063]

[0064] R i = p i × n i j i .

[0065] wherein r i is the coordinate vector of the observation satellite in the earth-centered inertial system, n i j is the direction vector obtained for the target observation, β is the included angle between the two vectors, is the estimated value of the position vector in the Kepler parameter, p i is the estimated value of the vector length from the observation satellite to the target satellite.

[0066] S104: according to the initial estimated value of the position vector obtained in S103, the target fitness function corresponding to the particle swarm algorithm is calculated, and the variables of the particle swarm algorithm are updated.

[0067] Further, the expression of the target fitness function in S104 is as follows:

[0068]

[0069] wherein R is the calculated position vector, F is the particle swarm fitness value corresponding to each iteration, and N is the number of observation times.

[0070] S105: the operations of S101 to S103 are repeated for iteration, the optimal target fitness function value is obtained, and the Kepler parameter initial value and the position vector value at each time are output.

[0071] Further, S2, the target position vector at each time obtained in S1 is orbit-determined, corrected, and outliers are eliminated by using the Laplace method and the least square algorithm, the orbit parameters corresponding to the short arc segment and the true anomaly corresponding to each time are calculated.

[0072] The formula used is as follows:

[0073] n = a -3 / 2 .

[0074] ΔM = n (t i -t0).

[0075]

[0076] p = a (1-e 2 ).

[0077] where M is the mean anomaly value, n is the average angular velocity of the orbit, a is the semi-major axis value of the orbit, e is the eccentricity, R is the position vector of the target, is the true anomaly.

[0078] And in order to obtain the eccentric anomaly change value ΔE, from the relationship between the mean anomaly and the eccentric anomaly Using Newton iteration method, the value of ΔE is obtained, and the iteration formula of E is as follows. After calculating the corresponding E value, difference is made, and the value of ΔE can be obtained:

[0079] Further, the specific steps of S2 are as follows, S201, using polynomial curve smoothing fitting, the position vector values at several time points are corrected.

[0080] Further, the specific content of the polynomial curve smoothing fitting in S201 is: obtaining the initial position vector and the velocity vector

[0081] Using polynomial fitting function for fitting.

[0082] The expression of the polynomial fitting function is: (j = 1, 2,..., n).

[0083] Wherein, R is the position vector value, t is the time value, d is the position vector at each time point, j is the position index,

[0084] S202: obtaining the position and velocity information at time t0, calculating the initial value of the semi-major axis, the mean anomaly value, the eccentric anomaly value, the parameter p and the parameter of the initial value.

[0085] S203: using the values obtained by S202, calculating the Laplace function value, and performing orbit extrapolation to obtain the position vector value and the velocity vector value at several time points.

[0086] Further, the expression of the Laplace function used in S203 is:

[0087] Wherein, V is the velocity of the target, f, g, F, G are Laplace operators.

[0088] Further, the orbit extrapolation expression in S203 is as follows:

[0089]

[0090] Wherein, a is the semi-major axis, r0 is the position vector module value, v0 is the velocity vector module value, p is the Laplace vector module value.

[0091] S204: Substitute the position vector values at several time instants into S201, repeat the contents of S201 to S203, and perform iterative calculation to obtain accurate values and accurate values, when and the iteration is terminated.

[0092] S205: Obtain several orbit roots using the position vector and velocity vector values and perform orbit determination.

[0093] The traditional short-arc orbit determination methods include the two-vector orbit determination algorithm and the Hibbert-Gibbs orbit determination algorithm.

[0094] The two-vector orbit determination algorithm uses the Gauss F-G function to solve the target orbit. When the measurement time interval is not large and the dynamic model is accurate, the velocity of the object is first roughly solved, and then the orbit parameters of the space object are calculated. The approximate algorithm is as follows. First, the approximate velocity of the space object is calculated

[0095]

[0096] The position vectors r1 and r2 measured at two time instants t1 and t2 are used to obtain v1 using the following formula.

[0097] v1 = [r2 - F(r1, v1, τ) r1] / G(r1, v1, τ).

[0098] Wherein:

[0099]

[0100] τ = t2 - t1.

[0101] Wherein F and G are Gauss functions of the orbit dynamics expansion.

[0102] Because the time interval is short, v1 is first used as the assumed value to calculate the F and G function values, and then a new v1 is calculated using the above formula, and new F and G function values are obtained by continuous calculation. The cycle is repeated until v1 meets the accuracy requirement. At this time, the orbit parameters of the object without maneuvering capability are calculated in the Earth inertial coordinate system. Although the two-vector orbit determination method uses less data, the accuracy of the orbit determination will be greatly deviated due to the inaccuracy of the model order, the deviation of the dynamic model, the disturbance, the observation error noise, etc.

[0103] The Hibbert-Gibbs orbit determination algorithm solves the initial position and velocity of the target through multiple vectors, and then solves the orbit parameters of the target. The specific formula implementation is as follows:

[0104]

[0105] wherein is the position, velocity vector value at each time.

[0106] After the specific position and velocity initial value is solved, the specific value of the orbit parameter is solved.

[0107] The following is a simulation of the short arc orbit determination algorithm. The target satellite is a near-earth satellite with a circular orbit, and the orbit radius is 7200 km. The observation frequency set for the target is 10 Hz, the observation time is 3 s, and there is a large observation error. The short arc orbit determination is performed on the target using the method, the double vector orbit determination and the Hibbert-Gibbs orbit determination algorithm. The results of the orbit determination are shown in the following table and Figure 2

[0108] Table 1 Comparison of orbit parameters obtained by each method

[0109]

[0110] Compared with the traditional orbit determination algorithm, the data range considered by the method is larger, and under a large observation error, the observation accuracy is maintained within an acceptable range, while the double vector orbit determination and the Hibbert-Gibbs orbit determination method have a large deviation due to the existence of the observation error, and the precision is greatly disturbed by the error. The method can effectively eliminate the influence of wild value, disturbance and the like. The number of sampling points and the sampling frequency can be increased to further improve the observation accuracy.

[0111] The orbit parameter value output by the short arc orbit determination can provide a better initial value for the precise orbit determination algorithm and orbit catalog identification.

[0112] Finally, it should be noted that the above examples are only used to illustrate the technical solutions of the present application and not to limit them. Although the present application has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present application can still be modified or replaced by equivalents, and these modifications or replacements cannot make the modified technical solutions deviate from the spirit and scope of the technical solutions of the present application.​

Claims

1. A Laplace shortest arc orbit determination algorithm based on particle swarm optimization, characterized in that: The following steps are involved: S1. Obtain the coordinate position vector of the observation satellite at each moment and its observed direction vector, and use the particle swarm algorithm to calculate the initial value of the Kepler parameter corresponding to the direction vector at each moment , and its corresponding target position vector in the inertial system, where is the semi-major axis value of the target orbit, Its eccentricity, is its near angle value; The specific steps of S1 are as follows: S101: Obtain three variables within the Kepler parameter range to obtain the particle swarm algorithm variables and define the three variables as the initial values ​​of the Kepler parameters ; S102: Calculating the mean near angle and its corresponding true near angle at each moment according to the Kepler parameter value, and calculating the modulus of the target position vector according to the mean near angle and its corresponding true near angle; S103: Obtain the observed direction vector and the coordinate position vector of the observed satellite and calculate the angle between the two to obtain an initial estimated value of the position vector of the target satellite; S104: Calculate the target fitness function corresponding to the particle swarm algorithm based on the initial estimated value of the position vector obtained in S103, and update the particle swarm algorithm variables; S105: Repeat the operations from S101 to S103, perform iterations, obtain the optimal target fitness function value, and output the initial value of the Kepler parameter and the position vector value at each moment; S2. Use the Laplace method and the least squares algorithm to determine the orbit, correct and eliminate outliers for the target position vector at each moment obtained in S1, and calculate the orbit parameters corresponding to the short arc segment and the true near angle corresponding to each moment. The specific steps of S2 are as follows: S201: using polynomial curve smoothing fitting to correct position vector values ​​at several moments; S202: Acquisition The position and velocity information at the moment is used to calculate the initial value of the orbital semi-major axis, the mean angle value, the eccentric angle value, and the parameters and parameters The initial value of S203: Using the values ​​obtained in S202, calculate the Laplace function value and perform orbit extrapolation to obtain the position vector value and velocity vector value at several moments; S204: Substitute the position vector values ​​at several moments into S201, repeat the contents of S201 to S203, perform iterative calculations, and obtain Exact value and Exact value; S205: Calculate a number of orbital elements using the position vector and velocity vector values ​​and perform orbit determination.

2. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 1, characterized in that: In step S102, the expression for calculating the near angle is as follows: ； in, is the near angle value, is the average angular velocity of the orbit; The calculation expression of the modulus is as follows: ； in, is the orbital semi-major axis value, is the eccentricity, is the target position vector, It is a position subscript.

3. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 2, characterized in that: The expression used to calculate the initial estimated value of the position vector of the target satellite in S103 is as follows: ； ； ； Among them, is the coordinate vector of the observation satellite in the geocentric inertial system, is the direction vector observed for the target, is the angle between the two vectors, is the estimated value of the position vector under Kepler parameters, is the estimated value of the radial modulus length from the observation satellite to the target satellite.

4. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 3, characterized in that: The expression of the objective fitness function in S104 is as follows: ； in, is the calculated position vector, is the particle swarm fitness value corresponding to each iteration, is the number of observation moments.

5. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 4, characterized in that: The specific content of polynomial curve smooth fitting in S201 is: Get the initial position vector at the moment , and the initial velocity vector ; Polynomial fitting function is used for fitting; The expression of the polynomial fitting function is: ； in, is the position vector value, is the time value, is the position vector at each moment, is the position subscript, .

6. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 1, characterized in that: The Laplace function expression used in calculating the Laplace function value in S203 is: ； in, is the target speed, is the Laplace operator.

7. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 6, characterized in that: The orbital extrapolation expression in S203 is: ； in, is the orbital semi-major axis value, is the modulus of the position vector, is the velocity vector modulus value, is the modulus of the Laplace vector.

8. The Laplace shortest arc orbit determination algorithm based on particle swarm optimization according to claim 7, characterized in that: During the iterative calculation process in S204, when as well as When , the iteration is terminated.

Images