CN114462256A - Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track - Google Patents

Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track Download PDF

Info

Publication number
CN114462256A
CN114462256A CN202210371828.2A CN202210371828A CN114462256A CN 114462256 A CN114462256 A CN 114462256A CN 202210371828 A CN202210371828 A CN 202210371828A CN 114462256 A CN114462256 A CN 114462256A
Authority
CN
China
Prior art keywords
vector
thrust
cooperative
moment
state quantity
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202210371828.2A
Other languages
Chinese (zh)
Other versions
CN114462256B (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.)
National University of Defense Technology
Original Assignee
National University of Defense Technology
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 National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN202210371828.2A priority Critical patent/CN114462256B/en
Publication of CN114462256A publication Critical patent/CN114462256A/en
Application granted granted Critical
Publication of CN114462256B publication Critical patent/CN114462256B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/242—Orbits and trajectories
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Chemical & Material Sciences (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Combustion & Propulsion (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

The application relates to a method, a device, computer equipment and a storage medium for determining a non-cooperative low-thrust maneuvering target track. The method comprises the following steps: the orbit determination method for the thrust acceleration modeling of the tracked non-cooperative target local orbit coordinate system is provided, an expansion state quantity transfer matrix calculation equation containing an acceleration component under a high-precision perturbation orbit model is deduced, orbit improvement is carried out on the non-cooperative small-thrust maneuvering target through radar observation data and a least square method, and finally orbit prediction is carried out through the determined expansion state quantity, so that the orbit determination of the non-cooperative target of continuous small-thrust maneuvering can be effectively carried out, and the method can be used for solving the problem of the orbit determination of the continuous small-thrust maneuvering of the non-cooperative satellite in a climbing or descending section.

Description

Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track
Technical Field
The application relates to the field of space situation perception, in particular to a method, a device, computer equipment and a storage medium for determining a non-cooperative low-thrust maneuvering target orbit of acceleration modeling under a local orbital system.
Background
The precise orbit determination and the high-precision orbit extrapolation of the non-cooperative maneuvering target have important significance on space situation perception tasks such as target cataloging, collision early warning and the like.
The star-chain satellite implements continuous small-thrust maneuvering, but the maneuvering strategy is unknown due to non-cooperation of the target, so that the thrust acceleration of the satellite is difficult to accurately model. The kinetic equation of the absolute orbital motion of the space target is generally described in a J2000 geocentric inertial system; however, in practical mission, in order to save energy, the orbital maneuver of the spacecraft is generally along the orbital trajectory or normal direction; because the spacecraft rotates around the earth and the direction vector of the tracking or normal maneuver in the inertial system also changes in a rotating way, the thrust acceleration vector is directly described in the inertial system, and a converged orbit improvement solution is difficult to obtain through state topology and least square estimation; therefore, the prior art is difficult to realize maneuvering section orbit determination and extrapolation prediction of the non-cooperative low-thrust maneuvering satellite.
Disclosure of Invention
In view of the above, there is a need to provide a method, an apparatus, a computer device and a storage medium for determining a non-cooperative low-thrust maneuvering target trajectory for acceleration modeling under a local track system, so as to improve the improvement and extrapolation effect of the low-thrust maneuvering target trajectory.
A method of non-cooperative low thrust maneuver target trajectory determination, the method comprising:
extracting station coordinates from observation data of a plurality of arc sections of a ground radar, and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment;
determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix;
according to the track state at the first moment, setting the initial acceleration value of the non-cooperative low-thrust maneuvering target at the first moment to be zero, and determining the initial expansion state quantity value of the non-cooperative low-thrust maneuvering target at the first moment;
obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula;
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model to predict the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
In one embodiment, the method further comprises the following steps: if the non-cooperative low-thrust maneuvering target has an initial value of cataloguing, acquiring the non-cooperative low-thrust maneuvering target from the cataloguing database
Figure 328887DEST_PATH_IMAGE001
Initial track state of time
Figure 495163DEST_PATH_IMAGE002
According to SGP4 method or other analytic perturbation orbit prediction model
Figure 367305DEST_PATH_IMAGE003
Obtaining the non-cooperative low-thrust maneuvering target
Figure 95089DEST_PATH_IMAGE004
Track state of time of day
Figure 583839DEST_PATH_IMAGE005
(ii) a Wherein the content of the first and second substances,
Figure 789693DEST_PATH_IMAGE006
is a position vector at a first time instant,
Figure 47499DEST_PATH_IMAGE007
is a velocity vector at a first time;
if the non-cooperative low-thrust maneuvering target has no cataloged initial value, a first observation point of a first radar observation arc section in the radar observation vector is observed
Figure 946185DEST_PATH_IMAGE008
And the last observation point
Figure 391072DEST_PATH_IMAGE009
Obtaining a non-cooperative maneuver object
Figure 400617DEST_PATH_IMAGE010
Position vector of time
Figure 778508DEST_PATH_IMAGE011
And
Figure 83981DEST_PATH_IMAGE012
position vector of time
Figure 281744DEST_PATH_IMAGE013
Further, furtherNon-cooperative maneuvering target obtained by adopting Lambert algorithm
Figure 94980DEST_PATH_IMAGE014
Track state of time of day
Figure 61799DEST_PATH_IMAGE015
(ii) a Wherein the subscript
Figure 302287DEST_PATH_IMAGE016
The total number of the observed data of the first observed arc segment.
In one embodiment, the method further comprises the following steps: the construction and description of the expansion state quantity of the non-cooperative small-thrust maneuvering target motion system are as follows:
Figure 987346DEST_PATH_IMAGE017
(ii) a Wherein the content of the first and second substances,
Figure 604272DEST_PATH_IMAGE018
respectively representing a position vector and a velocity vector under the geocentric inertial system;
Figure 425598DEST_PATH_IMAGE019
representing an acceleration vector under a local orbit coordinate system;
constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity, wherein the track dynamics model comprises the following steps:
Figure 102567DEST_PATH_IMAGE020
wherein the content of the first and second substances,
Figure 9343DEST_PATH_IMAGE021
is a constant of the gravity of the earth,
Figure 397337DEST_PATH_IMAGE022
the position velocity vector of the target under the geocentric inertial system is taken as a target;
Figure 604327DEST_PATH_IMAGE023
the acceleration of the non-spherical gravity perturbation of the earth,
Figure 921039DEST_PATH_IMAGE024
the third body gravity perturbation acceleration caused by the sun, the moon and other stars,
Figure 315111DEST_PATH_IMAGE025
in order to accelerate the air resistance,
Figure 273840DEST_PATH_IMAGE026
the solar light pressure perturbation acceleration is measured,
Figure 69758DEST_PATH_IMAGE027
representing perturbation acceleration due to tidal forces;
Figure 822950DEST_PATH_IMAGE028
for the thrust acceleration vector expressed in RTN,
Figure 969898DEST_PATH_IMAGE028
for modeling coefficients for acceleration based on current statistical models, e.g. it is advisable
Figure 201159DEST_PATH_IMAGE029
，
Figure 117162DEST_PATH_IMAGE030
For a one-step prediction of the acceleration, at each step of the integration,
Figure 536861DEST_PATH_IMAGE030
keeping the constant value of the input unchanged;
Figure 905526DEST_PATH_IMAGE031
is a transformation matrix from the local orbit coordinate system to the geocentric inertial system, which is expressed as:
Figure 206057DEST_PATH_IMAGE032
。
in one embodiment, the method further comprises the following steps: and performing partial derivative solution on the expansion state quantity by using a right function in the orbit dynamics model, and calculating a partial derivative matrix as follows:
Figure 976567DEST_PATH_IMAGE033
wherein the content of the first and second substances,
Figure 71562DEST_PATH_IMAGE034
a zero matrix of 3 rows and 3 columns is shown,
Figure 927522DEST_PATH_IMAGE035
an identity matrix, representing 3 rows and 3 columns
Figure 500586DEST_PATH_IMAGE036
、
Figure 125603DEST_PATH_IMAGE037
Is a matrix of partial derivatives of gravitational and perturbed acceleration components to a target position velocity vector,
Figure 657078DEST_PATH_IMAGE038
，
Figure 233291DEST_PATH_IMAGE039
matrix of
Figure 141204DEST_PATH_IMAGE040
、
Figure 620727DEST_PATH_IMAGE041
Partial derivative matrixes of target thrust acceleration components to target position velocity vectors respectively are all
Figure 323104DEST_PATH_IMAGE042
A matrix of (a);
determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix as follows:
Figure 888077DEST_PATH_IMAGE043
wherein, when the integral is solved by the above formula,
Figure 68523DEST_PATH_IMAGE044
the initial value is taken as
Figure 933711DEST_PATH_IMAGE045
Wherein
Figure 541409DEST_PATH_IMAGE046
Representing an identity matrix of 9 rows and 9 columns.
In one embodiment, the method further comprises the following steps: obtaining the extended state vector predicted values of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the extended state vector initial value at the first moment and the orbit dynamics model as follows:
Figure 859258DEST_PATH_IMAGE047
wherein the content of the first and second substances,
Figure 843395DEST_PATH_IMAGE048
representing an expansion state quantity forecast value of each radar observation data moment;
Figure 297510DEST_PATH_IMAGE049
the total number of the data observation time of the radar is;
obtaining a first-order state transition matrix of a plurality of radar observation data moments according to the initial value of the expanded state quantity at the first moment and the first-order state transition matrix formula:
Figure 341689DEST_PATH_IMAGE050
。
in one embodiment, the method further comprises the following steps: obtaining the radar observation vector as follows:
Figure 648299DEST_PATH_IMAGE051
wherein
Figure 170547DEST_PATH_IMAGE052
Respectively correspond to
Figure 744748DEST_PATH_IMAGE053
An observed value of a time;
obtaining observed quantity according to the radar observation vector
Figure 694250DEST_PATH_IMAGE054
And observation quantity to extended state quantity partial derivative matrix
Figure 721111DEST_PATH_IMAGE055
(ii) a Wherein the content of the first and second substances,
Figure 312630DEST_PATH_IMAGE056
，
Figure 475758DEST_PATH_IMAGE057
；
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value:
Figure 861740DEST_PATH_IMAGE058
calculating a partial derivative matrix from the state deviation to the observation residual according to the first-order state transition matrix as follows:
Figure 375898DEST_PATH_IMAGE059
。
and performing least square iteration improvement on the initial value of the expanded state quantity at the first moment according to the observation residual and a partial derivative matrix from the state deviation to the observation residual, wherein an iteration formula is as follows:
Figure 239949DEST_PATH_IMAGE060
wherein the content of the first and second substances,
Figure 523162DEST_PATH_IMAGE061
expressed in the second of least squares estimation iterationskStep (A) tok = 1,2,…,K) The target expansion state quantity of (2),Kis the preset maximum number of iterations,kwhen the ratio is not less than 1,
Figure 578581DEST_PATH_IMAGE062
(ii) a Let the azimuth standard difference observed by radar be
Figure 580035DEST_PATH_IMAGE063
Standard difference of pitch angle of
Figure 247777DEST_PATH_IMAGE064
Standard deviation of the skew distance of
Figure 385497DEST_PATH_IMAGE065
Then the weight matrix can be expressed as:
Figure 378861DEST_PATH_IMAGE066
when the iteration times are more than the preset times or the relative value of the root mean square error of the observed residual error is less than 10-6Then, finishing iteration and outputting the precision value of the expansion state quantity at the first moment
Figure 336452DEST_PATH_IMAGE067
。
In one embodiment, the method further comprises the following steps: the observations include azimuth, pitch, and skew data.
A non-cooperative low thrust maneuver target trajectory determination device, the device comprising:
the system comprises a radar observation vector acquisition module, a tracking module and a tracking module, wherein the radar observation vector acquisition module is used for extracting station coordinates from observation data of a plurality of arc sections of a foundation radar and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment; determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
the track dynamics model building module is used for building an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system and building a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
the first-order state transition matrix formula determining module is used for solving partial derivatives of the expanded state quantities by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantities according to the partial derivative matrix;
the least square iteration module is used for setting the initial acceleration value of the non-cooperative small-thrust maneuvering target at the first moment to be zero according to the track state at the first moment and determining the initial expansion state quantity value of the non-cooperative small-thrust maneuvering target at the first moment; obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula; calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
and the track prediction module is used for substituting the precision value of the expansion state quantity at the first moment into the track dynamics model to predict the flight track of the non-cooperative low-thrust maneuvering target at any moment.
A computer device comprising a memory and a processor, the memory storing a computer program, the processor implementing the following steps when executing the computer program:
extracting station coordinates from observation data of a plurality of arc sections of a ground radar, and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment;
determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix;
according to the track state at the first moment, setting the initial acceleration value of the non-cooperative low-thrust maneuvering target at the first moment to be zero, and determining the initial expansion state quantity value of the non-cooperative low-thrust maneuvering target at the first moment;
obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula;
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model to predict the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
A computer-readable storage medium, on which a computer program is stored which, when executed by a processor, carries out the steps of:
extracting station coordinates from observation data of a plurality of arc sections of a ground radar, and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment;
determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix;
according to the track state at the first moment, setting the initial acceleration value of the non-cooperative low-thrust maneuvering target at the first moment to be zero, and determining the initial expansion state quantity value of the non-cooperative low-thrust maneuvering target at the first moment;
obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula;
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model to predict the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
The method, the device, the computer equipment and the storage medium for determining the orbit of the non-cooperative low-thrust maneuvering target determine the radar observation vector of the tracked non-cooperative low-thrust maneuvering target from the observation data of a plurality of arc sections of the ground radar, and further determine the orbit state of the target at the first moment; constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix; obtaining an extended state vector prediction value of a non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value and the orbit dynamics model at the first moment, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value and a first-order state transition matrix formula at the first moment; performing track improvement on the non-cooperative low-thrust maneuvering target by radar observation data and a least square method until an iteration ending condition is met, and outputting an expansion state quantity precision value at a first moment; substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model, and predicting the flight trajectory of the non-cooperative small-thrust maneuvering target at any moment. The invention provides an orbit determination method for carrying out thrust acceleration modeling on a tracked non-cooperative target local orbit coordinate system, which deduces an extended state quantity transfer matrix calculation equation containing an acceleration component under a high-precision perturbation orbit model, carries out orbit improvement on a non-cooperative small-thrust maneuvering target through radar observation data and a least square method, and finally carries out orbit prediction through the determined extended state quantity, thereby effectively carrying out orbit determination on the continuous small-thrust maneuvering non-cooperative target.
Drawings
FIG. 1 is a schematic flow chart diagram of a method for determining a target trajectory for a non-cooperative low thrust maneuver in one embodiment;
FIG. 2 is a schematic flow chart diagram of a method for determining a target trajectory for a non-cooperative low thrust maneuver in an exemplary embodiment;
FIG. 3 is a block diagram of a non-cooperative low thrust maneuver target trajectory determination device in one embodiment;
FIG. 4 is a diagram illustrating an internal structure of a computer device according to an embodiment.
Detailed Description
In order to make the objects, technical solutions and advantages of the present application more apparent, the present application is described in further detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the present application and are not intended to limit the present application.
In one embodiment, as shown in fig. 1, there is provided a non-cooperative low thrust maneuver target trajectory determination method, comprising the steps of:
and 102, extracting station coordinates from observation data of a plurality of arc sections of the ground radar, and determining a radar observation vector of the tracked non-cooperative low-thrust maneuvering target according to the station coordinates.
Specifically, the site coordinates, the azimuth angle, the pitch angle and the slant range information of the observation non-cooperative maneuvering target are extracted from the observation data of a plurality of arc sections of the ground radar, and the site coordinates of the survey station S are recorded as [ lat, lon, alt ]]Wherein lat is the geographical latitude of the site, lon is the geographical longitude of the site, alt is the geographical elevation of the site; let us rememberiThe radar observation data is
Figure 73464DEST_PATH_IMAGE068
Wherein
Figure 331270DEST_PATH_IMAGE069
Is as followsiThe azimuth angle observed by the individual radar,
Figure 964377DEST_PATH_IMAGE070
is as followsiThe pitch angle as observed by the individual radar,
Figure 940423DEST_PATH_IMAGE071
is as followsiThe slant distance observed by each radar is obtained as the total observation vector
Figure 215547DEST_PATH_IMAGE072
Wherein
Figure 327859DEST_PATH_IMAGE073
Respectively correspond to
Figure 893052DEST_PATH_IMAGE074
The observed value of the time of day,Nthe total number of observed data points.
And 104, determining the track state of the non-cooperative low-thrust maneuvering target at the first moment according to the inventory database or the radar observation vector.
The orbital state includes a position vector and a velocity vector of the non-cooperative low-thrust maneuver object at the first time.
If the tracked non-cooperative maneuvering target has an initial value of the catalog, the non-cooperative maneuvering target is obtained according to the catalog database
Figure 90815DEST_PATH_IMAGE075
Initial track state of time
Figure 904050DEST_PATH_IMAGE076
Wherein
Figure 136448DEST_PATH_IMAGE077
For non-cooperative maneuvering objects
Figure 376937DEST_PATH_IMAGE078
The position vector of the time of day,
Figure 61996DEST_PATH_IMAGE079
for non-cooperative maneuvering objects
Figure 678922DEST_PATH_IMAGE075
A velocity vector of a time of day; the method adopts an SGP4 method or other analytic perturbation orbit forecasting models
Figure 500248DEST_PATH_IMAGE080
Forecasting the first radar observation data point moment
Figure 911637DEST_PATH_IMAGE081
Obtaining a non-cooperative maneuver at
Figure 818414DEST_PATH_IMAGE082
Track state of time of day
Figure 239031DEST_PATH_IMAGE083
。
If the tracked non-cooperative maneuvering target has no catalogued initial value, a first observation point based on a first radar observation arc segment
Figure 413398DEST_PATH_IMAGE084
And the last observation point
Figure 261268DEST_PATH_IMAGE085
Wherein the subscript
Figure 655340DEST_PATH_IMAGE086
Is the total number of observation data of the first observation arc segment, an
Figure 348490DEST_PATH_IMAGE087
Can solve the non-cooperative maneuvering target in
Figure 409987DEST_PATH_IMAGE088
Position vector of time
Figure 163179DEST_PATH_IMAGE089
Position vector of time
Figure 44547DEST_PATH_IMAGE090
Further adopting Lambert algorithm well known in the field of space dynamics to obtain non-cooperative maneuvering target
Figure 806967DEST_PATH_IMAGE091
Track state of time of day
Figure 457391DEST_PATH_IMAGE092
。
And 106, constructing an expansion state quantity describing the motion system of the non-cooperative small-thrust maneuvering target, and constructing a track dynamics model of the non-cooperative small-thrust maneuvering target according to the expansion state quantity.
The expansion state quantity comprises a position vector and a velocity vector under the geocentric inertial system and an acceleration vector under the local orbit coordinate system.
Since the maneuvering strategy and maneuvering acceleration of a non-cooperative low-thrust maneuvering target are unknown, three components of the acceleration vector under the target local orbit coordinate system are described as:
Figure 381485DEST_PATH_IMAGE093
wherein, in the step (A),
Figure 517193DEST_PATH_IMAGE094
、
Figure 552145DEST_PATH_IMAGE095
、
Figure 322655DEST_PATH_IMAGE096
respectively thrust acceleration in a target local orbit coordinate systemxA shaft,yA shaft,zThe component of the axis. The target local orbit coordinate system (hereinafter referred to as RTN system) is defined as the origin at the tracked target centroid,xthe axis is along the connecting line direction of the geocentric and the target barycenter,zthe shaft is normal to the track surface,ythe shaft and the other two shafts form a right-hand system. Further expanding the state quantity of the system to be estimated into
Figure 417650DEST_PATH_IMAGE097
In which
Figure 273611DEST_PATH_IMAGE098
The position vector of the target in the centroid inertia system (hereinafter referred to as ECI system),
Figure 112254DEST_PATH_IMAGE099
for a velocity vector targeted under the ECI series,
Figure 2850DEST_PATH_IMAGE100
the acceleration vector of the target under the RTN system is shown.
Modeling the orbit dynamics equation of the non-cooperative low-thrust maneuvering target as follows:
Figure 268746DEST_PATH_IMAGE101
Figure 612002DEST_PATH_IMAGE102
wherein the content of the first and second substances,
Figure 988757DEST_PATH_IMAGE103
is a constant of the gravity of the earth,
Figure 733859DEST_PATH_IMAGE104
position velocity vector under target ECI system;
Figure 669192DEST_PATH_IMAGE105
the acceleration of the non-spherical gravity perturbation of the earth,
Figure 234165DEST_PATH_IMAGE106
the third body gravity perturbation acceleration caused by the sun, the moon and other stars,
Figure 680190DEST_PATH_IMAGE107
in order to accelerate the air resistance,
Figure 279799DEST_PATH_IMAGE108
the solar light pressure perturbation acceleration is measured,
Figure 153077DEST_PATH_IMAGE109
representing perturbation acceleration due to tidal forces;
Figure 939767DEST_PATH_IMAGE110
for the thrust acceleration vector expressed in RTN,
Figure 189483DEST_PATH_IMAGE110
for modeling coefficients for acceleration based on current statistical models, e.g. it is advisable
Figure 909177DEST_PATH_IMAGE111
，
Figure 687778DEST_PATH_IMAGE112
For a one-step prediction of the acceleration, at each step of the integration,
Figure 227343DEST_PATH_IMAGE112
keeping the constant value of the input unchanged;
Figure 15171DEST_PATH_IMAGE113
the transformation matrix for RTN system to ECI system can be expressed as:
Figure 819398DEST_PATH_IMAGE114
Figure 34479DEST_PATH_IMAGE115
and 108, performing partial derivative solution on the expanded state quantity by using a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix.
Expanding state quantity of non-cooperative small-thrust maneuvering target
Figure 61340DEST_PATH_IMAGE116
The state transition equation is satisfied:
Figure 387280DEST_PATH_IMAGE117
wherein
Figure 81566DEST_PATH_IMAGE118
Is an initial
Figure 201969DEST_PATH_IMAGE119
The state of the time expansion is arbitrary in the futuretConstantly expanding a first-order state transition matrix of the state; obtaining initial position velocity according to step 1.2 or 1.3
Figure 716127DEST_PATH_IMAGE120
Then, can obtain
Figure 580178DEST_PATH_IMAGE121
Initial state of time
Figure 128971DEST_PATH_IMAGE122
If the iterative initial step acceleration is unknown, it can be taken
Figure 420275DEST_PATH_IMAGE123
(ii) a First order state transition matrix
Figure 421729DEST_PATH_IMAGE124
Satisfies the following conditions:
Figure 588006DEST_PATH_IMAGE125
Figure 991305DEST_PATH_IMAGE126
wherein the content of the first and second substances,
Figure 719090DEST_PATH_IMAGE127
the matrix being a differential equation of dynamics
Figure 942261DEST_PATH_IMAGE128
Middle right function pair expansion state quantity
Figure 679272DEST_PATH_IMAGE129
Partial derivative matrix, integral solution of ordinary differential equation
Figure 671499DEST_PATH_IMAGE130
Is taken as
Figure 570185DEST_PATH_IMAGE131
Wherein
Figure 546231DEST_PATH_IMAGE132
Representing an identity matrix of 9 rows and 9 columns.
Computing
Figure 555776DEST_PATH_IMAGE133
A matrix, which can be expressed as:
Figure 668088DEST_PATH_IMAGE134
Figure 737675DEST_PATH_IMAGE135
wherein the content of the first and second substances,
Figure 436903DEST_PATH_IMAGE136
a zero matrix of 3 rows and 3 columns is shown,
Figure 250139DEST_PATH_IMAGE137
representing an identity matrix of 3 rows and 3 columns,
Figure 216958DEST_PATH_IMAGE138
、
Figure 457446DEST_PATH_IMAGE139
can be calculated according to professional teaching materials or references related to the space dynamics, such as the following references: NASA Goddard Space Flight center. General Analysis Tool (GMAT) chemical specificities [ R]Greenbelt, MD 20771, month 6 2020.
Matrix array
Figure 142505DEST_PATH_IMAGE140
、
Figure 759431DEST_PATH_IMAGE141
Partial derivative matrixes of target thrust acceleration components and target position velocity vectors respectively, which are all
Figure 580757DEST_PATH_IMAGE142
Of the matrix of (a). By using superscriptsabRepresents the matrix ofaGo to the firstbElements of columns, e.g.
Figure 257726DEST_PATH_IMAGE143
Representation matrix
Figure 430081DEST_PATH_IMAGE144
Row 2, column 3 elements; note the book
Figure 850698DEST_PATH_IMAGE145
As components of the target position vector in the ECI system, i.e.
Figure 526530DEST_PATH_IMAGE146
(ii) a Note the book
Figure 374400DEST_PATH_IMAGE147
For the component of the target velocity vector in the ECI system, i.e.
Figure 267008DEST_PATH_IMAGE148
(ii) a Note book
Figure 225737DEST_PATH_IMAGE149
For the component of the target acceleration vector in the target RTN system, i.e.
Figure 287234DEST_PATH_IMAGE150
(ii) a Then the
Figure 774847DEST_PATH_IMAGE151
、
Figure 921794DEST_PATH_IMAGE152
The computational expression for each element is:
Figure 418635DEST_PATH_IMAGE153
Figure 334638DEST_PATH_IMAGE154
Figure 258732DEST_PATH_IMAGE155
Figure 627396DEST_PATH_IMAGE156
Figure 662348DEST_PATH_IMAGE157
Figure 698437DEST_PATH_IMAGE158
Figure 312475DEST_PATH_IMAGE159
Figure 168436DEST_PATH_IMAGE160
Figure 7079DEST_PATH_IMAGE161
Figure 897675DEST_PATH_IMAGE162
Figure 163571DEST_PATH_IMAGE163
Figure 506827DEST_PATH_IMAGE164
Figure 883582DEST_PATH_IMAGE165
Figure 628684DEST_PATH_IMAGE166
Figure 65482DEST_PATH_IMAGE167
Figure 896035DEST_PATH_IMAGE168
Figure 342059DEST_PATH_IMAGE169
Figure 174624DEST_PATH_IMAGE170
wherein intermediate variables are used
Figure 47902DEST_PATH_IMAGE171
The expression is as follows:
Figure 365751DEST_PATH_IMAGE172
Figure 349887DEST_PATH_IMAGE173
and step 110, according to the track state at the first moment, setting the initial acceleration value of the non-cooperative low-thrust maneuvering target at the first moment to be zero, and determining the initial expansion state quantity value of the non-cooperative low-thrust maneuvering target at the first moment.
At the first radar observation point
Figure 804002DEST_PATH_IMAGE174
At the moment, the initial acceleration value is set to zero
Figure 582603DEST_PATH_IMAGE175
。
And 112, obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state vector value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix at the plurality of radar observation data moments according to the initial extended state vector value at the first moment and a first-order state transition matrix formula.
According to the initial value
Figure 387748DEST_PATH_IMAGE176
At the same timeIntegral formula
Figure 175575DEST_PATH_IMAGE177
And formula
Figure 484197DEST_PATH_IMAGE178
Obtaining the expansion state of the target at each radar observation data moment
Figure 699277DEST_PATH_IMAGE179
And a first order state transition matrix at each radar observation time
Figure 991718DEST_PATH_IMAGE180
。
And step 114, calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from the state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on the initial value of the extended state quantity at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precise value at the first moment.
Extending state quantities by forecasted spacecraft
Figure 553543DEST_PATH_IMAGE181
Calculating an observation residual according to an observation equation
Figure 982250DEST_PATH_IMAGE182
And a matrix of partial derivatives of the state deviation to the observed residual
Figure 368232DEST_PATH_IMAGE183
。
Figure 882390DEST_PATH_IMAGE184
Figure 12020DEST_PATH_IMAGE185
Wherein, the radar observed quantity is calculated according to the state quantity under ECI system
Figure 295234DEST_PATH_IMAGE186
And observation quantity to extended state quantity partial derivative matrix
Figure 852117DEST_PATH_IMAGE187
The method can be obtained from aerospace dynamics professional textbooks or documents, and the invention is not repeated. For example, reference may be made to: theory and application of Liulin, Tangjing, satellite orbit]Beijing: electronics industry publishers, 2 nd edition 2015.
Expanding state quantities for non-cooperative low-thrust maneuvering targets
Figure 587992DEST_PATH_IMAGE188
Performing least square iteration improvement, wherein the iteration formula is as follows:
Figure 521313DEST_PATH_IMAGE189
Figure 659033DEST_PATH_IMAGE190
wherein, the azimuth standard difference observed by the radar is
Figure 386818DEST_PATH_IMAGE191
Standard difference of pitch angle of
Figure 875568DEST_PATH_IMAGE192
Standard deviation of the skew distance of
Figure 579957DEST_PATH_IMAGE193
Then the weight matrix can be expressed as:
Figure 837763DEST_PATH_IMAGE194
Figure 736449DEST_PATH_IMAGE195
and repeating the iteration until the iteration converges. For example, can be set as
Figure 446916DEST_PATH_IMAGE196
Or the relative value of the root mean square error of the observed residual is less than 10-6Time of flight
Figure 456460DEST_PATH_IMAGE197
Quitting, in which case the precise track state of the non-cooperative maneuvering target is obtained as
Figure 568772DEST_PATH_IMAGE198
(ii) a Otherwise, the least square iteration is unsuccessful, the method is not applicable any more, and a user needs to detect whether the tracked target is a continuous low-thrust maneuvering target or not, or whether the observed value is correct or not, or whether the initial value of the track iteration is reasonable or not.
And step 116, substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model, and predicting the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
Track condition to be obtained
Figure 638360DEST_PATH_IMAGE199
Substituting equation as initial value
Figure 836123DEST_PATH_IMAGE200
Integral prediction to arbitrarytAnd the flight track of the non-cooperative low-thrust maneuvering target in the observation blind area can be obtained at any moment, the movement direction of the target can be predicted according to the flight track, and possible threats of the target can be evaluated and early warned.
In particular, estimating states based on least squares
Figure 649358DEST_PATH_IMAGE201
The track extrapolation prediction is only suitable for the situation that the maneuvering strategy of the non-cooperative small-thrust maneuvering target is not changed, and if the target uses a thrust direction or magnitude different from that of the observation tracking section in the track extrapolation section, the track extrapolation accuracy is highWill be affected to different extents.
In the method for determining the orbit of the non-cooperative low-thrust maneuvering target, a radar observation vector of a tracked non-cooperative low-thrust maneuvering target is determined from observation data of a plurality of arc sections of a ground radar, and then the orbit state of the target at a first moment is determined; constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix; obtaining an extended state vector prediction value of a non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value and the orbit dynamics model at the first moment, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value and a first-order state transition matrix formula at the first moment; performing track improvement on the non-cooperative low-thrust maneuvering target by radar observation data and a least square method until an iteration ending condition is met, and outputting an expansion state quantity precision value at a first moment; substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model, and predicting the flight trajectory of the non-cooperative small-thrust maneuvering target at any moment. The invention provides an orbit determination method for performing thrust acceleration modeling on a tracked non-cooperative target local orbit coordinate system, which deduces an extended state quantity transfer matrix calculation equation containing an acceleration component under a high-precision perturbation orbit model, performs orbit improvement on a non-cooperative small-thrust maneuvering target through radar observation data and a least square method, and performs orbit prediction through the determined extended state quantity, so that the orbit determination can be effectively performed on the non-cooperative target of continuous small-thrust maneuvering, and the method can be used for solving the problem of performing the orbit determination of the continuous small-thrust maneuvering on a climbing or descending section of a non-cooperative satellite.
In one embodiment, the method further comprises the following steps: the observations include azimuth, pitch, and skew data.
It should be understood that, although the steps in the flowchart of fig. 1 are shown in order as indicated by the arrows, the steps are not necessarily performed in order as indicated by the arrows. The steps are not limited to being performed in the exact order illustrated and, unless explicitly stated herein, may be performed in other orders. Moreover, at least a portion of the steps in fig. 1 may include multiple sub-steps or multiple stages that are not necessarily performed at the same time, but may be performed at different times, and the order of performance of the sub-steps or stages is not necessarily sequential, but may be performed in turn or alternately with other steps or at least a portion of the sub-steps or stages of other steps.
In one embodiment, as shown in fig. 2, there is provided a method for determining a non-cooperative low thrust maneuver target trajectory, comprising:
step 1: and acquiring radar observation data and an initial orbit of the non-cooperative low-thrust maneuvering target.
1.1 extracting station address coordinates and information of azimuth angle, pitch angle and slant distance of observation non-cooperative maneuvering target from observation data of a plurality of arc sections of ground radar, wherein the station address coordinates of a survey station S are [46.8 deg,130.32 deg,101 m](ii) a Let us rememberiAn observation of radar is
Figure 881756DEST_PATH_IMAGE202
Wherein
Figure 617850DEST_PATH_IMAGE203
Is as followsiThe azimuth angle observed by the individual radar,
Figure 37330DEST_PATH_IMAGE204
is as followsiThe pitch angle as observed by the individual radar,
Figure 654256DEST_PATH_IMAGE205
is as followsiThe slant distance observed by each radar is obtained as the total observation vector
Figure 475582DEST_PATH_IMAGE206
Wherein
Figure 152551DEST_PATH_IMAGE207
Respectively correspond to
Figure 59327DEST_PATH_IMAGE208
The observed value of the time of day,Nthe total number of observed data points. And (4) switching to step 1.2 because the non-cooperative maneuvering target has an initial value of cataloguing.
1.2 obtaining non-cooperative maneuver targets from inventory database
Figure 948786DEST_PATH_IMAGE209
The initial orbit state at that time is:
Figure 155776DEST_PATH_IMAGE210
the method adopts an SGP4 method or other analytic perturbation orbit forecasting models
Figure 738067DEST_PATH_IMAGE211
Forecasting the first radar observation data point moment
Figure 132139DEST_PATH_IMAGE212
Obtaining a non-cooperative maneuver at
Figure 323824DEST_PATH_IMAGE213
The track state at the moment is:
Figure 119742DEST_PATH_IMAGE214
step 2: and establishing a high-precision perturbation orbit dynamics model of the non-cooperative small-thrust maneuvering target.
2.1 describes the three components of the acceleration vector in the target local orbit coordinate system as:
Figure 872934DEST_PATH_IMAGE215
wherein, in the step (A),
Figure 754302DEST_PATH_IMAGE216
respectively thrust acceleration in a target local orbit coordinate systemxA shaft,yA shaft,zComponent of the axis, the initial value of acceleration being taken as
Figure 516722DEST_PATH_IMAGE217
。
2.2 modeling the orbit dynamics equation of the non-cooperative low thrust maneuvering target as:
Figure 167146DEST_PATH_IMAGE218
Figure 91240DEST_PATH_IMAGE219
wherein the content of the first and second substances,
Figure 459904DEST_PATH_IMAGE220
is a constant of the gravity of the earth,
Figure 760436DEST_PATH_IMAGE221
for the acceleration modeling coefficient based on the current statistical model, take
Figure 530945DEST_PATH_IMAGE222
，
Figure 861826DEST_PATH_IMAGE223
For a one-step prediction of the acceleration, at each step of the integration,
Figure 983366DEST_PATH_IMAGE223
keeping the constant value of the input unchanged;
Figure 87588DEST_PATH_IMAGE224
the transformation matrix for RTN system to ECI system can be expressed as:
Figure 447025DEST_PATH_IMAGE225
Figure 978501DEST_PATH_IMAGE226
and step 3: and calculating a first-order state transition matrix of the expanded state quantity under the high-precision perturbation model.
3.1 expansion State quantity of non-cooperative Low thrust maneuver targets
Figure 56178DEST_PATH_IMAGE227
The state transition equation is satisfied:
Figure 698512DEST_PATH_IMAGE228
wherein
Figure 443614DEST_PATH_IMAGE229
Is an initial
Figure 880412DEST_PATH_IMAGE230
The state is expanded to any future at any momenttConstantly expanding a first-order state transition matrix of the state; obtaining initial position velocity according to step 1.2 or 1.3
Figure 710965DEST_PATH_IMAGE231
Then, can obtain
Figure 389945DEST_PATH_IMAGE232
Initial state of time
Figure 723975DEST_PATH_IMAGE233
If the iterative initial step acceleration is unknown, it can be taken
Figure 597253DEST_PATH_IMAGE234
(ii) a First order state transition matrix
Figure 649522DEST_PATH_IMAGE235
Satisfies the following conditions:
Figure 633659DEST_PATH_IMAGE236
Figure 353353DEST_PATH_IMAGE237
wherein the content of the first and second substances,
Figure 131953DEST_PATH_IMAGE238
the matrix being a differential equation of dynamics
Figure 937098DEST_PATH_IMAGE239
Middle right function pair expansion state quantity
Figure 724926DEST_PATH_IMAGE240
Partial derivative matrix, integral solution of ordinary differential equation
Figure 33547DEST_PATH_IMAGE241
Is taken as
Figure 248628DEST_PATH_IMAGE242
Wherein
Figure 275490DEST_PATH_IMAGE243
An identity matrix of 9 rows and 9 columns is represented,
Figure 97034DEST_PATH_IMAGE244
the matrix is calculated by step 3.2.
3.2 calculation of
Figure 260163DEST_PATH_IMAGE244
A matrix, which can be expressed as:
Figure 646145DEST_PATH_IMAGE245
Figure 425882DEST_PATH_IMAGE246
wherein the content of the first and second substances,
Figure 555512DEST_PATH_IMAGE247
a zero matrix of 3 rows and 3 columns is shown,
Figure 573146DEST_PATH_IMAGE248
representing an identity matrix of 3 rows and 3 columns,
Figure 395609DEST_PATH_IMAGE249
can be calculated according to professional teaching materials or references related to the space dynamics, such as the following references: NASA Goddard Space Flight center. General Analysis Tool (GMAT) chemical specificities [ R]Greenbelt, MD 20771, month 6 2020.
3.3 matrix
Figure 397063DEST_PATH_IMAGE250
Partial derivative matrixes of target thrust acceleration components and target position velocity vectors respectively, which are all
Figure 64805DEST_PATH_IMAGE251
Of the matrix of (a).
And 4, step 4: orbit determination is performed based on a weighted least squares estimate.
4.1 at the first observation Point
Figure 202525DEST_PATH_IMAGE252
At the moment, the initial value of the acceleration is
Figure 195889DEST_PATH_IMAGE253
Target initial position velocity obtained according to step 1.2 or 1.3
Figure 917595DEST_PATH_IMAGE254
In the second of least squares estimation iterationkStep (A) tok = 1,2,…,K) Recording the target expansion state as
Figure 389027DEST_PATH_IMAGE255
WhereinKIs the set maximum number of iterations.kWhen the ratio is not less than 1,
Figure 646833DEST_PATH_IMAGE256
。
4.2 according to the initial value
Figure 545519DEST_PATH_IMAGE257
Simultaneous integral formula
Figure 255986DEST_PATH_IMAGE258
And formula
Figure 265531DEST_PATH_IMAGE259
Obtaining the expansion state of the target at each radar observation data moment
Figure 643422DEST_PATH_IMAGE260
And a first order state transition matrix at each radar observation time
Figure 713009DEST_PATH_IMAGE261
。
4.3 expanding State quantities by forecasted spacecraft
Figure 910773DEST_PATH_IMAGE262
Calculating an observation residual according to an observation equation
Figure 724008DEST_PATH_IMAGE263
And a matrix of partial derivatives of the state deviation to the observed residual
Figure 690827DEST_PATH_IMAGE264
。
Figure 931315DEST_PATH_IMAGE265
Figure 117839DEST_PATH_IMAGE266
Wherein, the radar observed quantity is calculated according to the state quantity under ECI system
Figure 734766DEST_PATH_IMAGE267
And observation quantity to extended state quantity partial derivative matrix
Figure 556091DEST_PATH_IMAGE268
The method can be obtained from aerospace dynamics professional textbooks or documents, and the invention is not repeated. For example, reference may be made to: theory and application of Liulin, Tangjing, satellite orbit]Beijing: electronics industry publishers, 2 nd edition 2015.
4.4 expansion State quantities for non-cooperative Low-thrust maneuver targets
Figure 233060DEST_PATH_IMAGE269
Performing least square iteration improvement, wherein the iteration formula is as follows:
Figure 405415DEST_PATH_IMAGE270
wherein the azimuth standard difference observed by the radar is
Figure 294874DEST_PATH_IMAGE271
Standard difference of pitch angle of
Figure 501864DEST_PATH_IMAGE272
Standard deviation of slope distance of
Figure 84155DEST_PATH_IMAGE273
Then the weight matrix can be expressed as:
Figure 478228DEST_PATH_IMAGE274
4.5 repeat steps 4.2 to 4.4 until the iteration converges. For example, can be set as
Figure 436956DEST_PATH_IMAGE275
Or the relative value of the root mean square error of the observed residual error is less than 1e-6
Figure 232874DEST_PATH_IMAGE276
Quitting, wherein the initial value of the precise track of the non-cooperative maneuvering target is obtained
Figure 750181DEST_PATH_IMAGE277
(ii) a Otherwise, the least square iteration is unsuccessful, the method is not applicable any more, and a user needs to detect whether the tracked target is a continuous low-thrust maneuvering target or not, or whether the observed value is correct or not, or whether the initial value of the track iteration is reasonable or not.
And 5: orbit state based on least square estimation
Figure 631549DEST_PATH_IMAGE278
And extrapolating and forecasting the future flight track of the non-cooperative low-thrust maneuvering target.
5.1 track status obtained in step 4
Figure 128390DEST_PATH_IMAGE279
Substituting equation as initial value
Figure 44393DEST_PATH_IMAGE280
Integral prediction to arbitrarytAnd the flight track of the non-cooperative low-thrust maneuvering target in the observation blind area can be obtained at any moment, the movement direction of the target can be predicted according to the flight track, and possible threats of the target can be evaluated and early warned.
5.2 in particular, estimating the state based on least squares
Figure 968487DEST_PATH_IMAGE281
The track extrapolation prediction is only applicable to the situation that the maneuvering strategy of the non-cooperative small-thrust maneuvering target is not changed, and if the target uses a thrust direction or magnitude different from that of the observation tracking section in the track extrapolation section, the track extrapolation accuracy is influenced to different degrees.
Determined by the above steps
Figure 602730DEST_PATH_IMAGE282
Forecasting the track state to the moment of the last observation point to obtain the track state as follows:
Figure 637682DEST_PATH_IMAGE283
using the value as an initial value, substituting the initial value into an equation
Figure 408192DEST_PATH_IMAGE284
And (4) performing orbit prediction, performing 12 click-through from 54968.566 seconds to 30 days in 6 months in 2021, and comparing with the real TLE orbit data of the satellite with the star link number 48465 at the moment to obtain prediction accuracy. Obtaining forecast track as
Figure 768767DEST_PATH_IMAGE285
By looking up the TLE orbit data of the satellite with the number of 48465 satellite chains, the real orbit of the satellite at the same moment is obtained
Figure 624727DEST_PATH_IMAGE286
Alternatively, the satellite orbit may be predicted based on a conventional free orbit extrapolation method, and the prediction errors calculated for the method of the present invention and the conventional free extrapolation method are shown in table 1.
TABLE 1 prediction accuracy of the method of the present invention after determination of the orbit of 48465 satellite chain satellites
Figure 463370DEST_PATH_IMAGE287
As can be seen from Table 1, the method can effectively determine the orbit of the non-cooperative maneuvering target, and the forecasting precision of the determined orbit more than half a day is obviously superior to that of the traditional method without considering the influence of the orbit maneuvering.
Therefore, the application example can prove that the method for determining the orbit of the non-cooperative small-thrust maneuvering target, provided by the invention, can effectively model and represent the acceleration influence of the orbit determination problem of the non-cooperative maneuvering target, and has higher orbit determination precision.
In one embodiment, as shown in fig. 3, there is provided a non-cooperative low-thrust maneuvering target trajectory determination device, including: a radar observation vector obtaining module 302, a trajectory dynamics model establishing module 304, a first order state transition matrix formula determining module 306, a least squares iteration module 308, and a trajectory prediction module 310, wherein:
a radar observation vector obtaining module 302, configured to extract site coordinates from observation data of multiple arc segments of a ground-based radar, and determine a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the site coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative low-thrust maneuvering target at all times; determining the track state of the non-cooperative low-thrust maneuvering target at the first moment according to the cataloging database or the radar observation vector; the orbit state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuvering target at a first moment;
the track dynamics model establishing module 304 is used for establishing an expansion state quantity for describing a motion system of the non-cooperative low-thrust maneuvering target and establishing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
a first-order state transition matrix formula determining module 306, configured to perform partial derivative solution on the expanded state quantity by using a right function in the orbit dynamics model to obtain a partial derivative matrix, and determine a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix;
the least square iteration module 308 is configured to determine an initial value of an expansion state quantity of the non-cooperative low-thrust maneuvering target at the first moment according to the track state at the first moment and by setting an initial value of acceleration of the non-cooperative low-thrust maneuvering target at the first moment to zero; obtaining an extended state vector prediction value of a non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value and the orbit dynamics model at the first moment, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value and a first-order state transition matrix formula at the first moment; calculating an observation residual error according to a radar observation vector and an expanded state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to a first-order state transition matrix, performing least square iteration improvement on an expanded state quantity initial value at a first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an expanded state quantity precise value at the first moment;
and the trajectory prediction module 310 is configured to substitute the precision value of the expansion state quantity at the first moment into the orbit dynamics model, and predict the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
The radar observation vector acquisition module 302 is further configured to acquire the non-cooperative low-thrust maneuvering target from the inventory database if the non-cooperative low-thrust maneuvering target has an inventory initial value
Figure 353966DEST_PATH_IMAGE288
Initial track state of time
Figure 115467DEST_PATH_IMAGE289
According to SGP4 method or other analytic perturbation orbit prediction model
Figure 458724DEST_PATH_IMAGE290
Obtain a non-cooperative low-thrust maneuvering target
Figure 101058DEST_PATH_IMAGE291
Track state of time of day
Figure 846160DEST_PATH_IMAGE292
(ii) a Wherein the content of the first and second substances,
Figure 282958DEST_PATH_IMAGE293
is a position vector at a first time instant,
Figure 113510DEST_PATH_IMAGE294
is a velocity vector at a first time; if the non-cooperative low-thrust maneuvering target has no cataloged initial value, observing a first observation point of a first radar observation arc section in a radar observation vector
Figure 28377DEST_PATH_IMAGE295
And finallyOne observation point
Figure 893565DEST_PATH_IMAGE296
Obtaining a non-cooperative maneuver object
Figure 766843DEST_PATH_IMAGE297
Position vector of time
Figure 819112DEST_PATH_IMAGE298
And
Figure 803249DEST_PATH_IMAGE299
position vector of time
Figure 21478DEST_PATH_IMAGE300
Further adopting Lambert algorithm to obtain non-cooperative maneuvering target
Figure 800079DEST_PATH_IMAGE297
Track state of time of day
Figure 605223DEST_PATH_IMAGE301
(ii) a Wherein the subscript
Figure 127472DEST_PATH_IMAGE302
The total number of the observed data of the first observed arc segment.
The orbit dynamics model building module 304 is further configured to build an expansion state quantity describing a non-cooperative small-thrust maneuvering target motion system as follows:
Figure 701672DEST_PATH_IMAGE303
(ii) a Wherein the content of the first and second substances,
Figure 651174DEST_PATH_IMAGE304
respectively representing a position vector and a velocity vector under the geocentric inertial system;
Figure 678036DEST_PATH_IMAGE305
representing an acceleration vector under a local orbit coordinate system;
the track dynamics model of the non-cooperative low-thrust maneuvering target is constructed according to the expansion state quantity as follows:
Figure 269554DEST_PATH_IMAGE306
。
wherein the content of the first and second substances,
Figure 432682DEST_PATH_IMAGE307
is a constant of the gravity of the earth,
Figure 818664DEST_PATH_IMAGE308
the position velocity vector of the target under the geocentric inertial system is taken as a target;
Figure 332822DEST_PATH_IMAGE309
the acceleration of the non-spherical gravity perturbation of the earth,
Figure 698338DEST_PATH_IMAGE310
the third body gravity perturbation acceleration caused by the sun, the moon and other stars,
Figure 981552DEST_PATH_IMAGE311
in order to accelerate the air resistance,
Figure 804014DEST_PATH_IMAGE312
the solar light pressure perturbation acceleration is measured,
Figure 539889DEST_PATH_IMAGE313
representing perturbation acceleration due to tidal forces;
Figure 207631DEST_PATH_IMAGE314
for the thrust acceleration vector expressed in RTN,
Figure 610930DEST_PATH_IMAGE314
for modeling coefficients for acceleration based on current statistical models, e.g. it is advisable
Figure 338715DEST_PATH_IMAGE315
，
Figure 561886DEST_PATH_IMAGE316
For a one-step prediction of the acceleration, at each step of the integration,
Figure 33318DEST_PATH_IMAGE316
keeping the constant value of the input unchanged;
Figure 291124DEST_PATH_IMAGE317
is a transformation matrix from the local orbit coordinate system to the geocentric inertial system, which is expressed as:
Figure 189810DEST_PATH_IMAGE318
。
the first-order state transition matrix formula determining module 306 is further configured to perform partial derivative solution on the expanded state quantities by using a right function in the orbit dynamics model, and calculate a partial derivative matrix as:
Figure 398812DEST_PATH_IMAGE319
wherein the content of the first and second substances,
Figure 673936DEST_PATH_IMAGE320
a zero matrix of 3 rows and 3 columns is shown,
Figure 520669DEST_PATH_IMAGE321
an identity matrix, representing 3 rows and 3 columns
Figure 855836DEST_PATH_IMAGE322
、
Figure 53599DEST_PATH_IMAGE323
Is a matrix of partial derivatives of gravitational and perturbed acceleration components to a target position velocity vector,
Figure 866834DEST_PATH_IMAGE324
，
Figure 99232DEST_PATH_IMAGE325
matrix of
Figure 74141DEST_PATH_IMAGE326
Partial derivative matrixes of target thrust acceleration components to target position velocity vectors respectively are all
Figure 24780DEST_PATH_IMAGE327
A matrix of (a);
the first-order state transition matrix formula for determining the expanded state quantity according to the partial derivative matrix is as follows:
Figure 641706DEST_PATH_IMAGE328
wherein, when the integral is solved by the above formula,
Figure 463031DEST_PATH_IMAGE329
the initial value is taken as
Figure 874421DEST_PATH_IMAGE330
Wherein
Figure 112740DEST_PATH_IMAGE331
Representing an identity matrix of 9 rows and 9 columns.
The least square iteration module 308 is further configured to obtain an extended state vector predicted value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state vector value at the first moment and the orbit dynamics model as follows:
Figure 533357DEST_PATH_IMAGE332
wherein the content of the first and second substances,
Figure 209189DEST_PATH_IMAGE333
representing an expansion state quantity forecast value of each radar observation data moment;
Figure 57060DEST_PATH_IMAGE334
the total number of times of observing data for the radar;
according to the initial value of the expanded state quantity at the first moment and a first-order state transition matrix formula, obtaining a first-order state transition matrix at a plurality of radar observation data moments as follows:
Figure 451132DEST_PATH_IMAGE335
。
the least squares iteration module 308 is further configured to obtain a radar observation vector as:
Figure 144281DEST_PATH_IMAGE336
wherein
Figure 205778DEST_PATH_IMAGE337
Respectively correspond to
Figure 958971DEST_PATH_IMAGE338
An observed value of a time;
obtaining observed quantity according to radar observation vector
Figure 574760DEST_PATH_IMAGE339
And observation quantity to extended state quantity partial derivative matrix
Figure 570135DEST_PATH_IMAGE340
(ii) a Wherein the content of the first and second substances,
Figure 486139DEST_PATH_IMAGE341
，
Figure 551178DEST_PATH_IMAGE342
；
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value:
Figure 123105DEST_PATH_IMAGE343
calculating a partial derivative matrix from the state deviation to the observation residual according to the first-order state transition matrix as follows:
Figure 128363DEST_PATH_IMAGE344
。
according to the observation residual error and a partial derivative matrix from the state deviation to the observation residual error, performing least square iteration improvement on the initial value of the expanded state quantity at the first moment, wherein the iteration formula is as follows:
Figure 633294DEST_PATH_IMAGE345
wherein the content of the first and second substances,
Figure 728289DEST_PATH_IMAGE346
representing iteration in least squares estimationkStep (A) tok = 1,2,…,K) The target expansion state quantity of (2),Kis the maximum iteration number which is preset,kwhen the ratio is not less than 1,
Figure 318670DEST_PATH_IMAGE347
(ii) a Let the azimuth standard difference observed by radar be
Figure 157313DEST_PATH_IMAGE348
Standard difference of pitch angle of
Figure 47909DEST_PATH_IMAGE349
Standard deviation of the skew distance of
Figure 579384DEST_PATH_IMAGE350
Then the weight matrix can be expressed as:
Figure 657062DEST_PATH_IMAGE351
when the iteration times are more than the preset times or the relative value of the root mean square error of the observed residual error is less than 1e-6Then, finishing iteration and outputting the precision value of the expansion state quantity at the first moment
Figure 797931DEST_PATH_IMAGE352
。
For specific limitations of the non-cooperative low-thrust maneuvering target trajectory determination device, reference may be made to the above limitations of the non-cooperative low-thrust maneuvering target trajectory determination method, and details thereof are not repeated here. The various modules in the non-cooperative low-thrust maneuvering target trajectory determination device can be realized in whole or in part by software, hardware and a combination thereof. The modules can be embedded in a hardware form or independent from a processor in the computer device, and can also be stored in a memory in the computer device in a software form, so that the processor can call and execute operations corresponding to the modules.
In one embodiment, a computer device is provided, which may be a terminal, and its internal structure diagram may be as shown in fig. 4. The computer device includes a processor, a memory, a network interface, a display screen, and an input device connected by a system bus. Wherein the processor of the computer device is configured to provide computing and control capabilities. The memory of the computer device comprises a nonvolatile storage medium and an internal memory. The non-volatile storage medium stores an operating system and a computer program. The internal memory provides an environment for the operation of an operating system and computer programs in the non-volatile storage medium. The network interface of the computer device is used for communicating with an external terminal through a network connection. The computer program is executed by a processor to implement a non-cooperative low thrust maneuver target trajectory determination method. The display screen of the computer equipment can be a liquid crystal display screen or an electronic ink display screen, and the input device of the computer equipment can be a touch layer covered on the display screen, a key, a track ball or a touch pad arranged on the shell of the computer equipment, an external keyboard, a touch pad or a mouse and the like.
Those skilled in the art will appreciate that the architecture shown in fig. 4 is merely a block diagram of some of the structures associated with the disclosed aspects and is not intended to limit the computing devices to which the disclosed aspects apply, as particular computing devices may include more or less components than those shown, or may combine certain components, or have a different arrangement of components.
In an embodiment, a computer device is provided, comprising a memory storing a computer program and a processor implementing the steps of the above method embodiments when executing the computer program.
In an embodiment, a computer-readable storage medium is provided, on which a computer program is stored, which computer program, when being executed by a processor, carries out the steps of the above-mentioned method embodiments.
It will be understood by those skilled in the art that all or part of the processes of the methods of the embodiments described above may be implemented by hardware instructions of a computer program, which may be stored in a non-volatile computer-readable storage medium, and when executed, may include the processes of the embodiments of the methods described above. Any reference to memory, storage, database, or other medium used in the embodiments provided herein may include non-volatile and/or volatile memory, among others. Non-volatile memory can include read-only memory (ROM), Programmable ROM (PROM), Electrically Programmable ROM (EPROM), Electrically Erasable Programmable ROM (EEPROM), or flash memory. Volatile memory can include Random Access Memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in a variety of forms such as Static RAM (SRAM), Dynamic RAM (DRAM), Synchronous DRAM (SDRAM), Double Data Rate SDRAM (DDRSDRAM), Enhanced SDRAM (ESDRAM), Synchronous Link DRAM (SLDRAM), Rambus Direct RAM (RDRAM), direct bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM).
The technical features of the above embodiments can be arbitrarily combined, and for the sake of brevity, all possible combinations of the technical features in the above embodiments are not described, but should be considered as the scope of the present specification as long as there is no contradiction between the combinations of the technical features.
The above-mentioned embodiments only express several embodiments of the present application, and the description thereof is more specific and detailed, but not construed as limiting the scope of the invention. It should be noted that, for a person skilled in the art, several variations and modifications can be made without departing from the concept of the present application, which falls within the scope of protection of the present application. Therefore, the protection scope of the present patent shall be subject to the appended claims.

Claims (10)

1. A method for determining a non-cooperative low thrust maneuver target trajectory, the method comprising:
extracting station coordinates from observation data of a plurality of arc sections of a ground radar, and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment;
determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
constructing an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system, and constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
performing partial derivative solving on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix;
according to the track state at the first moment, setting the initial acceleration value of the non-cooperative low-thrust maneuvering target at the first moment to be zero, and determining the initial expansion state quantity value of the non-cooperative low-thrust maneuvering target at the first moment;
obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula;
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
substituting the precision value of the expansion state quantity at the first moment into the orbit dynamics model to predict the flight trajectory of the non-cooperative low-thrust maneuvering target at any moment.
2. The method of claim 1, wherein the orbital state of the non-cooperative low-thrust maneuver target at the first time is determined from a cataloged database or the radar observation vectors; the trajectory state includes a position vector and a velocity vector of the non-cooperative low-thrust maneuver object at a first time, including:
if the non-cooperative low-thrust maneuvering target has an initial value of cataloguing, acquiring the non-cooperative low-thrust maneuvering target from the cataloguing database
Figure 558758DEST_PATH_IMAGE001
Initial track state of time
Figure 944740DEST_PATH_IMAGE002
According to SGP4 method or other analytic perturbation orbit prediction model
Figure 193319DEST_PATH_IMAGE003
Obtaining the non-cooperative low-thrust maneuvering target
Figure 322949DEST_PATH_IMAGE004
Track state of time of day
Figure 871742DEST_PATH_IMAGE005
(ii) a Wherein the content of the first and second substances,
Figure 163046DEST_PATH_IMAGE006
is a position vector at a first time instant,
Figure 164500DEST_PATH_IMAGE007
is a velocity vector at a first time;
if the non-cooperative low-thrust maneuvering target has no cataloged initial value, a first observation point of a first radar observation arc section in the radar observation vector is observed
Figure 97821DEST_PATH_IMAGE008
And the last observation point
Figure 501120DEST_PATH_IMAGE009
Obtaining a non-cooperative maneuver object
Figure 482369DEST_PATH_IMAGE010
Position vector of time
Figure 971119DEST_PATH_IMAGE011
And
Figure 176972DEST_PATH_IMAGE012
position vector of time
Figure 434778DEST_PATH_IMAGE013
Further adopting Lambert algorithm to obtain non-cooperative maneuvering target
Figure 333464DEST_PATH_IMAGE014
Of time of dayTrack condition
Figure 309510DEST_PATH_IMAGE015
(ii) a Wherein the subscript
Figure 319055DEST_PATH_IMAGE016
The total number of the observed data of the first observed arc segment.
3. The method according to claim 2, characterized by constructing an expansion state quantity describing a non-cooperative small-thrust maneuvering target moving system, and constructing an orbit dynamics model of the non-cooperative small-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system, and comprises the following steps:
the construction and description of the expansion state quantity of the non-cooperative small-thrust maneuvering target motion system are as follows:
Figure 696946DEST_PATH_IMAGE017
(ii) a Wherein the content of the first and second substances,
Figure 500954DEST_PATH_IMAGE018
、
Figure 698717DEST_PATH_IMAGE019
respectively representing a position vector and a velocity vector under the geocentric inertial system;
Figure 511953DEST_PATH_IMAGE020
representing an acceleration vector under a local orbit coordinate system;
constructing a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity, wherein the track dynamics model comprises the following steps:
Figure 744351DEST_PATH_IMAGE021
wherein the content of the first and second substances,
Figure 483374DEST_PATH_IMAGE022
is a constant of the gravity of the earth,
Figure 168434DEST_PATH_IMAGE023
、
Figure 519781DEST_PATH_IMAGE024
the position velocity vector of the target under the geocentric inertial system is taken as a target;
Figure 341106DEST_PATH_IMAGE025
the acceleration of the non-spherical gravity perturbation of the earth,
Figure 752496DEST_PATH_IMAGE026
the third body gravity perturbation acceleration caused by the sun, the moon and other stars,
Figure 924851DEST_PATH_IMAGE027
in order to accelerate the air resistance,
Figure 814310DEST_PATH_IMAGE028
the solar light pressure perturbation acceleration is measured,
Figure 21300DEST_PATH_IMAGE029
representing perturbation acceleration due to tidal forces;
Figure 869171DEST_PATH_IMAGE030
for the thrust acceleration vector expressed in RTN,
Figure 263243DEST_PATH_IMAGE031
for the acceleration modeling coefficients based on the current statistical model,
Figure 723436DEST_PATH_IMAGE032
for one-step prediction of acceleration, at each step of integration，
Figure 253775DEST_PATH_IMAGE033
Keeping the constant value of the input unchanged;
Figure 272546DEST_PATH_IMAGE034
is a transformation matrix from the local orbit coordinate system to the geocentric inertial system, which is expressed as:
Figure 153915DEST_PATH_IMAGE035
。
4. the method according to claim 3, wherein the partial derivative solution is performed on the expanded state quantity by a right function in the orbit dynamics model to obtain a partial derivative matrix, and a first-order state transition matrix formula of the expanded state quantity is determined according to the partial derivative matrix, including:
and performing partial derivative solution on the expansion state quantity by using a right function in the orbit dynamics model, and calculating a partial derivative matrix as follows:
Figure 916334DEST_PATH_IMAGE036
wherein the content of the first and second substances,
Figure 832338DEST_PATH_IMAGE037
a zero matrix of 3 rows and 3 columns is shown,
Figure 490852DEST_PATH_IMAGE038
an identity matrix, representing 3 rows and 3 columns
Figure 125096DEST_PATH_IMAGE039
、
Figure 160048DEST_PATH_IMAGE040
Is a pair of gravitational and gravitational acceleration componentsA matrix of partial derivatives of the target position velocity vector,
Figure 196137DEST_PATH_IMAGE041
，
Figure 25553DEST_PATH_IMAGE042
matrix of
Figure 645628DEST_PATH_IMAGE043
、
Figure 749850DEST_PATH_IMAGE044
Partial derivative matrixes of target thrust acceleration components to target position velocity vectors respectively are all
Figure 374866DEST_PATH_IMAGE045
A matrix of (a);
determining a first-order state transition matrix formula of the expanded state quantity according to the partial derivative matrix as follows:
Figure 640763DEST_PATH_IMAGE046
wherein, when the integral is solved by the above formula,
Figure 718440DEST_PATH_IMAGE047
the initial value is taken as
Figure 626353DEST_PATH_IMAGE048
Wherein
Figure 371455DEST_PATH_IMAGE049
Representing an identity matrix of 9 rows and 9 columns.
5. The method of claim 4, wherein obtaining extended state vector prediction values of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix at the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula comprises:
obtaining the extended state vector predicted values of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the extended state vector initial value at the first moment and the orbit dynamics model as follows:
Figure 808253DEST_PATH_IMAGE050
wherein the content of the first and second substances,
Figure 638806DEST_PATH_IMAGE051
representing an expansion state quantity forecast value of each radar observation data moment;
Figure 84831DEST_PATH_IMAGE052
the total number of the data observation time of the radar is;
obtaining a first-order state transition matrix of a plurality of radar observation data moments according to the initial value of the expanded state quantity at the first moment and the first-order state transition matrix formula:
Figure 418860DEST_PATH_IMAGE053
。
6. the method of claim 5, wherein calculating an observation residual according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from a state deviation to the observation residual according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual and the partial derivative matrix from the state deviation to the observation residual until an iteration end condition is met, and outputting an extended state quantity precise value at the first moment comprises:
obtaining the radar observation vector as follows:
Figure 292138DEST_PATH_IMAGE054
wherein
Figure 840013DEST_PATH_IMAGE055
Respectively correspond to
Figure 558570DEST_PATH_IMAGE056
An observed value of a time;
obtaining observed quantity according to the radar observation vector
Figure 278265DEST_PATH_IMAGE057
And observation quantity to extended state quantity partial derivative matrix
Figure 322444DEST_PATH_IMAGE058
(ii) a Wherein the content of the first and second substances,
Figure 862010DEST_PATH_IMAGE059
，
Figure 649837DEST_PATH_IMAGE060
；
calculating an observation residual error according to the radar observation vector and the extended state vector predicted value:
Figure 224038DEST_PATH_IMAGE061
calculating a partial derivative matrix from the state deviation to the observation residual according to the first-order state transition matrix as follows:
Figure 439119DEST_PATH_IMAGE062
and performing least square iteration improvement on the initial value of the expanded state quantity at the first moment according to the observation residual and a partial derivative matrix from the state deviation to the observation residual, wherein an iteration formula is as follows:
Figure 200401DEST_PATH_IMAGE063
wherein the content of the first and second substances,
Figure 526340DEST_PATH_IMAGE064
representing iteration in least squares estimationkStep (A) tok = 1,2,…,K) The target expansion state quantity of (2),Kis the preset maximum number of iterations,kwhen the pressure is not greater than 1, the pressure is lower than 1,
Figure 955048DEST_PATH_IMAGE065
(ii) a Let the azimuth standard difference observed by radar be
Figure 839565DEST_PATH_IMAGE066
Standard difference of pitch angle of
Figure 353723DEST_PATH_IMAGE067
Standard deviation of the skew distance of
Figure 483353DEST_PATH_IMAGE068
Then the weight matrix can be expressed as:
Figure 766567DEST_PATH_IMAGE069
when the iteration times are more than the preset times or the relative value of the root mean square error of the observed residual error is less than 10-6Then, finishing iteration and outputting the precision value of the expansion state quantity at the first moment
Figure 323450DEST_PATH_IMAGE070
。
7. The method of any of claims 1-6, wherein the observations comprise azimuth, pitch, and roll data.
8. A non-cooperative low thrust maneuver target trajectory determination device, comprising:
the system comprises a radar observation vector acquisition module, a tracking module and a tracking module, wherein the radar observation vector acquisition module is used for extracting station coordinates from observation data of a plurality of arc sections of a foundation radar and determining a radar observation vector of a tracked non-cooperative low-thrust maneuvering target according to the station coordinates; the radar observation vector comprises observation values of a plurality of radar observation data to the non-cooperative small-thrust maneuvering target at the moment; determining the orbit state of the non-cooperative low-thrust maneuvering target at a first moment according to a cataloging database or the radar observation vector; the orbital state comprises a position vector and a velocity vector of the non-cooperative low-thrust maneuver target at a first time;
the track dynamics model building module is used for building an expansion state quantity for describing a non-cooperative low-thrust maneuvering target motion system and building a track dynamics model of the non-cooperative low-thrust maneuvering target according to the expansion state quantity; the expansion state quantity comprises a position vector and a velocity vector under a geocentric inertial system and an acceleration vector under a local orbit coordinate system;
the first-order state transition matrix formula determining module is used for solving partial derivatives of the expanded state quantities by a right function in the orbit dynamics model to obtain a partial derivative matrix, and determining a first-order state transition matrix formula of the expanded state quantities according to the partial derivative matrix;
the least square iteration module is used for determining an initial value of the expansion state quantity of the non-cooperative low-thrust maneuvering target at the first moment according to the track state at the first moment and setting the initial value of the acceleration of the non-cooperative low-thrust maneuvering target at the first moment to be zero; obtaining an extended state vector forecast value of the non-cooperative low-thrust maneuvering target at a plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the orbit dynamics model, and obtaining a first-order state transition matrix of the plurality of radar observation data moments according to the initial extended state quantity value at the first moment and the first-order state transition matrix formula; calculating an observation residual error according to the radar observation vector and the extended state vector predicted value, calculating a partial derivative matrix from state deviation to the observation residual error according to the first-order state transition matrix, performing least square iteration improvement on an initial extended state quantity value at the first moment according to the observation residual error and the partial derivative matrix from the state deviation to the observation residual error until an iteration ending condition is met, and outputting an extended state quantity precision value at the first moment;
and the track prediction module is used for substituting the precision value of the expansion state quantity at the first moment into the track dynamics model to predict the flight track of the non-cooperative low-thrust maneuvering target at any moment.
9. A computer device comprising a memory and a processor, the memory storing a computer program, wherein the processor implements the steps of the method of any one of claims 1 to 7 when executing the computer program.
10. A computer-readable storage medium, on which a computer program is stored, which, when being executed by a processor, carries out the steps of the method of any one of claims 1 to 7.
CN202210371828.2A 2022-04-11 2022-04-11 Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track Active CN114462256B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210371828.2A CN114462256B (en) 2022-04-11 2022-04-11 Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210371828.2A CN114462256B (en) 2022-04-11 2022-04-11 Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track

Publications (2)

Publication Number Publication Date
CN114462256A true CN114462256A (en) 2022-05-10
CN114462256B CN114462256B (en) 2022-06-21

Family

ID=81417725

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210371828.2A Active CN114462256B (en) 2022-04-11 2022-04-11 Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track

Country Status (1)

Country Link
CN (1) CN114462256B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116108319A (en) * 2023-04-10 2023-05-12 中国人民解放军32035部队 Orbit forecasting method for constant thrust mode continuous maneuvering satellite

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8174433B1 (en) * 2009-02-26 2012-05-08 Raytheon Company Bias estimation and orbit determination
EP3514571A1 (en) * 2018-01-18 2019-07-24 Thales Method for tracking an aerial target, and radar implementing such a method
CN110567461A (en) * 2019-08-01 2019-12-13 北京航空航天大学 Non-cooperative spacecraft attitude and parameter estimation method considering no gyroscope
CN111994304A (en) * 2020-08-31 2020-11-27 北京理工大学 Low-thrust long-term position keeping method for geostationary orbit satellite
CN112084581A (en) * 2020-09-24 2020-12-15 中国人民解放军国防科技大学 Spacecraft low-thrust perturbation intersection trajectory optimization method and system
US20210026524A1 (en) * 2019-07-25 2021-01-28 ExoAnalytic Solutions, Inc. Systems and visualization interfaces for orbital paths and path parameters of space objects
CN113031638A (en) * 2021-03-10 2021-06-25 哈尔滨工业大学 Constellation deployment method, device and storage medium based on time-sharing phase modulation of low-thrust satellite

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8174433B1 (en) * 2009-02-26 2012-05-08 Raytheon Company Bias estimation and orbit determination
EP3514571A1 (en) * 2018-01-18 2019-07-24 Thales Method for tracking an aerial target, and radar implementing such a method
US20210026524A1 (en) * 2019-07-25 2021-01-28 ExoAnalytic Solutions, Inc. Systems and visualization interfaces for orbital paths and path parameters of space objects
CN110567461A (en) * 2019-08-01 2019-12-13 北京航空航天大学 Non-cooperative spacecraft attitude and parameter estimation method considering no gyroscope
CN111994304A (en) * 2020-08-31 2020-11-27 北京理工大学 Low-thrust long-term position keeping method for geostationary orbit satellite
CN112084581A (en) * 2020-09-24 2020-12-15 中国人民解放军国防科技大学 Spacecraft low-thrust perturbation intersection trajectory optimization method and system
CN113031638A (en) * 2021-03-10 2021-06-25 哈尔滨工业大学 Constellation deployment method, device and storage medium based on time-sharing phase modulation of low-thrust satellite

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
ZHEN YANG 等: "Adaptive Finite-Time Control on SE(3) for Spacecraft Final Proximity Maneuvers with Input Quantization", 《• INTERNATIONAL JOURNAL OF AEROSPACE ENGINEERING》 *
刘磊: "空间目标轨道机动探测研究", 《中国优秀硕士学位论文全文数据库（自然地理学和测绘学）》 *
崔红正等: "不同推力下的非合作空间目标轨道机动检测", 《宇航学报》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116108319A (en) * 2023-04-10 2023-05-12 中国人民解放军32035部队 Orbit forecasting method for constant thrust mode continuous maneuvering satellite
CN116108319B (en) * 2023-04-10 2023-08-11 中国人民解放军32035部队 Orbit forecasting method for constant thrust mode continuous maneuvering satellite

Also Published As

Publication number Publication date
CN114462256B (en) 2022-06-21

Similar Documents

Publication Publication Date Title
Hu et al. Model predictive based unscented Kalman filter for hypersonic vehicle navigation with INS/GNSS integration
March et al. High-fidelity geometry models for improving the consistency of CHAMP, GRACE, GOCE and Swarm thermospheric density data sets
Hajiyev et al. Gyro-free attitude and rate estimation for a small satellite using SVD and EKF
Li et al. Radio beacons/IMU integrated navigation for Mars entry
Park et al. Estimating small-body gravity field from shape model and navigation data
Pérez et al. Neural Network based calibration of atmospheric density models
CN110378012B (en) Strict regression orbit design method, system and medium considering high-order gravity field
CN105737823A (en) GPS/SINS/CNS integrated navigation method based on five-order CKF
CN101059349A (en) Minitype combined navigation system and self-adaptive filtering method
CN107402903B (en) Nonlinear system state deviation evolution method based on differential algebra and Gaussian sum
CN106767837A (en) Based on the spacecraft attitude method of estimation that volume quaternary number is estimated
CN108875244B (en) Orbit prediction precision improvement method based on random forest
CN105224737A (en) A kind of extraterrestrial target improvement of orbit just value correction method
Li et al. Innovative Mars entry integrated navigation using modified multiple model adaptive estimation
CN110059439B (en) Spacecraft orbit determination method based on data driving
CN114462256B (en) Method, device, equipment and medium for determining non-cooperative low-thrust maneuvering target track
CN115523927B (en) GEO spacecraft maneuvering detection method based on optical sensor observation
Singh et al. Feasibility of quasi-frozen, near-polar and extremely low-altitude lunar orbits
Visser et al. Horizontal and vertical thermospheric cross-wind from GOCE linear and angular accelerations
Williams et al. Validation of Linear Covariance Techniques for Mars Entry, Descent, and Landing Guidance and Navigation Performance Analysis
CN115265532A (en) Auxiliary filtering method for marine integrated navigation
Zhu et al. Fast approximation of optimal perturbed long-duration impulsive transfers via artificial neural networks
Fang et al. State estimation of double-pyramid tethered satellite formations using only two GPS sensors
CN111814313A (en) Method for designing regression orbit in high-precision gravitational field
Ananthasayanam Tuning of the Kalman filter using constant gains

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