CN116363206B - Single-star multistage trajectory missile track reconstruction method based on straight line reconstruction method - Google Patents

Abstract

The invention discloses a single-star multistage trajectory missile track reconstruction method based on a straight line reconstruction method, which comprises the steps of converting target two-dimensional coordinates of infrared images of continuous time sequences and early-warning satellite three-dimensional coordinates into a unified coordinate system, and extracting a time sequence satellite observation vector; estimating the firing angle of the multi-stage ballistic missile according to a firing estimation algorithm; reconstructing a discrete three-dimensional track of a multi-stage ballistic missile target by adopting a linear reconstruction method according to the satellite observation vector of the continuous time sequence and the calculated shot; constructing acceleration characteristics and rationality limiting conditions of the multistage ballistic missile, and screening suspected trajectories; and carrying out track smoothing filtering on the discrete multi-stage ballistic missile track by using the volume Kalman filtering with multiple weight judgment to obtain a smooth multi-stage ballistic missile three-dimensional track point with the same length as the observation vector. The method can reconstruct the track of the multi-stage ballistic missile according to few prior information of the multi-stage ballistic missile.

Description

Single-star multistage trajectory missile track reconstruction method based on straight line reconstruction method

Technical Field

The invention relates to the technical fields of photogrammetry, linear reconstruction and three-dimensional track reconstruction, in particular to a single-star multistage trajectory missile track reconstruction method based on a linear reconstruction method.

Background

Because the engine can be started only when the ballistic missile target is in the active section and observed by the infrared camera, the infrared sensor for observing the satellite, which is configured by the early warning system, only passively detects, namely only acquires the two angle values of the sight of the ballistic missile target, and cannot acquire depth signals. In order to solve the difficulty, students at home and abroad conduct a plurality of researches, and complement the imperfection caused by a single star observation positioning mode according to various prior information. Current research modes based on prior information are divided into two categories: the method is characterized in that modeling is carried out according to the motion characteristics of the missile, then a tracking task is completed according to the model, and the method is called a track feature modeling method, and the method is called a trajectory track template method by recording track priori feature information of various missiles to construct an information template library.

However, the multiple iterations of the track feature modeling method do not solve the problem that the track reconstruction accuracy gradually decreases with the increase of the reconstruction time due to various error factors, and even diverges finally. The method for acquiring the trajectory template is to pay more attention to acquiring more trajectory missile information or to ingeniously inductive analyze different trajectory missile trajectory characteristics, and the information of the missiles belongs to information which is difficult to acquire, so that the trajectory template library is easy to have insufficient information of the trajectory of the missiles, and serious errors are generated when the missiles which are not recorded are observed.

Therefore, there is a need to solve the above-mentioned problems.

Disclosure of Invention

The invention aims to: the invention aims to provide a single-star multistage trajectory missile trajectory reconstruction method based on a linear reconstruction method, which can calculate a three-dimensional trajectory of a target according to a time sequence space-based image by the linear reconstruction method, and finally uses a volume particle filtering algorithm to carry out trajectory smoothing.

The technical scheme is as follows: in order to achieve the above purpose, the invention discloses a single-star multistage trajectory missile track reconstruction method based on a straight line reconstruction method, which comprises the following steps:

(1) Converting the target two-dimensional coordinates of the infrared images of the continuous time sequence and the three-dimensional coordinates of the early warning satellite into a unified coordinate system, and extracting a time sequence satellite observation vector;

(2) According to the obtained observation vector, estimating the firing angle of the multi-stage ballistic missile according to a firing estimation algorithm;

(3) Reconstructing a discrete three-dimensional track of a multi-stage ballistic missile target by adopting a linear reconstruction method according to the satellite observation vector of the continuous time sequence and the calculated shot;

(4) Constructing acceleration characteristics and rationality limiting conditions of the multi-stage ballistic missile, and screening suspected trajectories according to the acceleration characteristics and rationality limiting conditions of the ballistic missile;

(5) And carrying out track smoothing filtering on the discrete multi-stage ballistic missile track by using the volume Kalman filtering with multiple weight judgment to obtain a smooth multi-stage ballistic missile three-dimensional track point with the same length as the observation vector.

Wherein step (1) comprises the steps of:

(1.1) converting a target two-dimensional coordinate based on a pixel coordinate system in a satellite image into an image plane coordinate system according to a conversion mode of the pixel coordinate system and the image plane coordinate system, wherein the conversion mode of the two is as follows:

, wherein ,/>Is the position of the target in the image plane coordinate system; />Is the position of the target in the pixel coordinate system; />Is the width of a two-dimensional photo,/>Is the line height of a two-dimensional photograph,/>The size of the pixels of the two-dimensional photo can be obtained from the camera parameters of the satellite;

(1.2) converting the target coordinates in the image plane coordinate system into the satellite body coordinate system according to the linear projection conversion of the image plane coordinate system and the satellite body coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>Is the object ofA position in a satellite body coordinate system; />Representative is the focal length of the sensor, which can be obtained from the camera parameters of the satellite;

(1.3) converting the target coordinates in the satellite body coordinate system into the satellite orbit measurement coordinate system according to the conversion matrix of the satellite body coordinate system and the satellite orbit measurement coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>For the position of the object in the orbit measurement coordinate system,for the position of the object under the body coordinate system, +.>Is based on the attitude angle->The expression mode of the created rotation matrix is as follows:

, wherein ,/>For observing the pitch angle, the roll angle and the yaw angle of the satellite, the pitch angle, the roll angle and the yaw angle of the satellite can be obtained according to a measurement device built in the satellite;

(1.4) converting the target coordinates in the satellite orbit measurement coordinate system into the WGS84 fixed coordinate system according to a conversion matrix of the satellite orbit measurement coordinate system and the WGS84 fixed coordinate system, wherein the conversion mode is as follows:

, wherein ,/>For satellites in the WGS84 geodetic fixed coordinate system,for the target in the WGS84 fixed coordinate system,/L>Is according to the latitude of the satellite point below>With longitude->The expression mode of the created rotation matrix is as follows:

(1.5) set->For the launch time point of a multi-stage ballistic missile target, < + >>Establishing an observation vector set for a theoretical shutdown point of a target object of the multi-stage ballistic missile:wherein is a->Time->Obtained by the formula:

, wherein />Is a proportional coefficient->The orientation element in the camera can be obtained from satellite camera parameters; the observation vector at any time is obtained from the relational expression.

Preferably, step (2) comprises the steps of:

(2.1) coordinates of three-dimensional image points in the WGS84 coordinate systemAnd observation vector->The constructed straight line and radius are +.>Spherical intersection of sphere radius +.>For the earth radius->Estimated altitude with ballistic missile target>And (3) summing; intersection of straight line and sphere->For the three-dimensional coordinate point corresponding to the target emission point, < >>The elements of the formula are as follows:

(2.2) calculating according to the formula of the step (2.1) to obtain the intersection pointAbout sphere radius->The functional expression is->The missile target corresponding to the observation times is expressed as +.>Obtaining the sphere radius of multiple observation moments according to the formula of the step (2.1)>The functional relation of (2) is +.>；

(2.3) parameters due to the orientation of the planeVery small approximately 0, & gt>And (3) withTwo points solve a plane of projection +.>And is directed to plane->The normal vector of (2) is:

；

(2.4) orderAnd->The defined plane of incidence is +.>Then->And->The plane of the shot is defined asThe method comprises the steps of carrying out a first treatment on the surface of the Altering +.in the formula of step (2.3)>2 and->When, i.e. calculate the acquisition plane +.>Normal vector of->And plane->Normal vector of->Plane +.>And plane->Normal vector of->And->Point at the same point and/or>And->Coincide and all pass through the sphere center of the earth and +.>Point, origin of geocentric coordinates +. >The following relation is given:

，

in the formula Is->About->、/>Functional representation of (a) and the sameIs->About->、/>A functional representation of (a);

(2.5) taking three points in the observation process，/>，/>I.e. in the formula of step (2.4)>The value of (2) is 3, and the formula of the step (2.4) is solved by using the Monte Carlo method to obtain +.>，/>，/>Solution of (2); the coordinates of the multi-stage ballistic missile launching point under the WGS84 coordinate system are +.>The normal vector of the coil plane passing through the missile launching point is:

，

in the formula Is equal to +.about.in WGS84 coordinate system>The axes are in the same direction and the direction vector modulo 1, then the shot is:

。

furthermore, the step (3) includes the steps of:

(3.1) constructing a camera projection matrix for projecting a multi-stage ballistic missile onto a photographCamera projection matrix->Has the following definition formula:

, wherein />For pixel coordinate system->Scale factor of axis>For pixel coordinate system->Scale factor of axis>As a scale factor, < >>Is the coordinates of the center of the photo in the pixel coordinate system, < >>Coordinates projected in pixel coordinate system for two-dimensional object,/->Is the coordinates of the object in the geodetic fixed coordinate system,

for the coordinate rotation matrix, +.>For coordinate translation matrix, +.>A projection matrix called a camera;

(3.2) assume that along a straight linePoint projection of a mobile multi-stage ballistic missile target +. >And let->Is straight line->In view->Projection onto, due to->Straight line in the image plane +.>On the other hand, the flow of qi is->At this time, a straight line +.>Let view->The expression in the WGS84 solid coordinate system is:

wherein parameter->、/>、/>View +.>Time observation line of sight +.>、/>、；/>Can be substituted into camera coordinates, namely satellite coordinates>Obtaining;

(3.3) setting a straight lineThe parametric equation expression of (2) is:

, wherein />、/>、/>For view->The three-dimensional coordinates of the moving target at the corresponding moment;

along a certain straight line to the targetThree-dimensional position ∈>And straight line->Direction vector->The following constraint formula is set:

，

，

wherein Normal vector to the warp plane of the target emission point, +.>Is constant (I)>Refers to a priori velocity information of the target; so far, all three-dimensional straight line spaces are embedded in a five-dimensional projection space and are constrained by the constraint formula;

(3.4) passing through a straight lineAnd perpendicular to view->Is>The method comprises the following steps:

thus->The expression of (2) can be solved by:

after rotation of the coordinate system, the straight line is represented under the pixel coordinate system>：

Its parameterization is expressed as +>Is provided with->Is straight line->Go up to any point and letFor->In view->Projection onto; obviously- >Because of->Straight line in the image plane +.>Go up and->；

According to camera projection matrix at different momentsAnd a time-series three-dimensional coordinate point of the object->For unknown straight line->The following linear constraint formula is established:

, wherein />For straight line along the target->Projection on photo, parameterized as +.>The method comprises the steps of carrying out a first treatment on the surface of the When->The unique solution can be solved at this time, and when +.>When the characteristic difference between the actual track and the straight line is too large or more than five views can obtain a least square solution; in general, the least squares estimation matrix has a rank of 5 and +.>To estimate the zero space of the matrix; when->When the track of the camera projection center is straight during certain tracking periods of the camera motion, the situation that the estimated matrix rank is 4 can occur, and the constraint is acted at the moment; at this time there is a ∈>、/>Extended two-dimensional zero space, thus->Scalar->The optimal value can be calculated according to the constraint formula in the step (3.3);

(3.5) according to straight linesRecursion of parameter equation, reconstructing discrete coordinates of the moving point, and solving straight line according to the linear constraint formula in the step (3.4)>Expressed by the following parametric equation:

，

wherein Is straight line->Go up to a point, let go of>Is straight line->Direction parameters of (2);

Straight lineIs>The solution can be found by:

；

(3.6) recording all suspected track results:

, wherein />Is a suspected trace->Three-dimensional coordinates of the time instant>Is->Time of moment->Taking the value according to the three-dimensional coordinate number>And (5) taking a value according to the number of suspected trajectories.

Further, the step (4) includes the steps of:

(4.1) establishing an engine and mass change model, wherein the engine thrust model is as follows:

, wherein ,/>Is the amount of time, representing the moment of movement of the target; />Is the elevation of the target at the current moment; />Is the end time of different active segment stages; />The fuel consumption rate, which is the target of the different active stage phases, and in most cases +.>Is a definite value, can be regarded as->；/>Is the stable exhaust speed of the gas nozzle at different active stage, and can be considered as +.>Is a fixed value; />Is the cross sectional area of the engine nozzle at different active stage; />Is the static pressure of air flow at the nozzles of different active stage; />Is the atmospheric static pressure of the elevation where the target flying process is located; the ballistic missiles of the same model are not very different>、/>、/>、/>、/>The values are only needed to be modified according to a transmitting program in actual early warning, so that the time sequence thrust information of a certain missile can be obtained;

The mass change model is as follows:

, wherein ,/>Is the initial mass of different active segment stages, and is generally the missile with the same model，/>、/>The value of (2) is fixed;

(4.2) modeling non-engine acceleration, air resistanceThe general expression of (2) is:

, wherein />Is dynamic pressure (is->Is air density->Is the velocity of the missile relative to air; />Is the maximum cross-sectional area of the missile body; />For Mach number>Record->Is the speed of sound; />As the drag coefficient, the footnote "0" indicates the attack angle of the missile +.>Is 0;

due to the air resistance of the missile target and the air resistance of the missile targetThe included angle of the movement direction is a flat angle, so that the resistance accelerationOpposite to the direction vector of the motion direction, can be expressed as follows:

, wherein ,/>Numerical representation of the characteristics of the missile shape design>As a function of the air density as a function of the elevation, wherein +.>，；/>Is a resistance function, and can represent the force born by the missile according to the speed value;

namely, the analysis and calculation of acceleration caused by the gravity of the earth to the ballistic missile target are realized, the numerical value of the gravity is calculated by using the expression of the space sphere, and the potential energy function obtained by the numerical value is taken asThe method comprises the steps of carrying out a first treatment on the surface of the Let the earth's long half shaft be +. >The short half shaft is +.>And the mass of the earth is symmetrical about the equatorial plane and all longitudinal planes, then the gravitational function of the earth can be expressed as:

, wherein />For the distance between the earth center and the target->For the quality of the object, +.>Gravity inclination angle->For perturbation term is constant, +.>Stress constant->Is Legendre polynomial;

in practice perturbation itemsThe values are not large, so in the earth gravity calculation, get +.>The perturbation term, the gravity function expression at this time, is:

, wherein ,/>；

Order theCorresponding gravitational acceleration->Radial in the earth>North gravityInclination angle->Is as follows:

projecting gravitational acceleration to the geocentric sagittal approach>Rotation axis of earthThe following steps are:

；

(4.3) setting an active segment ballistic base constraint,

and (4.4) setting specialized constraint aiming at the ballistic missile target, and carrying out suspected track screening according to the rationality constraint of the ballistic missile target consisting of the initiative section ballistic foundation constraint and the specialized constraint aiming at the ballistic missile target, wherein the acceleration constraint and the rationality constraint together, so that the error track section can be eliminated.

Preferably, the active segment ballistic basic constraints in step (4.3) include:

outside conditions: the trajectory of a ballistic missile cannot be within the earth, i.e. the elevation of each point on the trajectory needs to be greater than 0, assuming For the elevation of each point on the missile track, the constraint is as follows:

，

elevation increment: the elevation of points on the trajectory of the ballistic missile engine during the working phase must be consistent with the characteristic of increasing with time, and the specific constraint expression is as follows:

，

speed extremum: the maximum scalar rate of points on the missile target trajectory cannot exceed the first cosmic speed because the missile target cannot reach such a high speed with the power of the missile engine, and the first cosmic speed is recorded asThe actual constraint expression is as follows:

，

rate change with time: the scalar rate at each moment on the trajectory of the ballistic missile engine working phase is required to strictly conform to the gradually increasing constraint, and the actual constraint expression is as follows:

acceleration extremum: the maximum scalar acceleration at each moment on the trajectory of the ballistic missile engine working phase cannot be greater than +.>The specific constraint expression is as follows:

，

the initial point elevation condition is the initial detection point elevation condition: after azimuth information is obtained according to the observation sight of the initial point, longitude and latitude near the initial point can be calculated, and then elevation of the initial point can be calculated according to the local cloud layer elevation information at the momentConstraint range of->；

Gas injection limit: the fuel gas injection speed determines the thrust force of the missile engine, and two different forms of fuel injection rate constraints are set, wherein the liquid fuel injection rate is 2500-4600m/s, and the solid injection rate is 2000-3000m/s;

The viewing constraint: the distance between the observation vector and the intersection point of the earth cannot be further than the distance between the missile target and the satellite.

Furthermore, the specialized constraints for the ballistic missile target in step (4.4) include:

active segment vertical ascent segment durationConstraint:

，

the active section rises vertically Duan GaochengConstraint:

，

first-stage boosting transverse displacementConstraint:

，

duration of the whole active sectionConstraint:

，

whole active section flying elevationConstraint:

，

flight transverse-longitudinal ratio of whole driving sectionConstraint:

。

further, step (5) comprises the steps of:

(5.1) the volume Kalman filtering starts to select sigma points according to the estimated state, finally obtains some sigma points with the same weight, and then carries out state estimation and estimation covariance calculation according to the set system model prediction and iteration state values;

and (5.2) performing volume particle filtering, namely obtaining state prediction at each moment by adopting a state estimation and covariance estimation calculation process in the volume Kalman filtering step (5.1), taking a covariance matrix of the state prediction as a mean value and a variance in a particle filtering algorithm to establish a proposal distribution function, and obtaining a smooth multistage ballistic missile three-dimensional track point with the same length as an observation vector after volume particle filtering iterative calculation.

Preferably, step (5.1) comprises the steps of:

(5.1.1) initializing state and covariance, setting vector length of state equation，/>Linear reconstructed trajectory data for tracking a preliminary stage ballistic missile target>Calculating the state and covariance matrix of the missile:

，/>, wherein />To solve the expectations;

(5.1.2) time update, calculationPersonal->The point:

，/>, in the formula ：/>Is->Personal->Point (S)>Is thatDecomposing;

(5.1.3) propagationPersonal->Point, one-step prediction of state:

, wherein />In discrete form of a state equation;

inputting a matrix for a process;

and is also provided with，/>Constant determined for the target a priori acceleration means, < +.>Time intervals for satellite observations;

at this time, one-step prediction of the state is regarded as an average value of the state at this time:

, wherein />Expressed as pair->、/>、/>Coordinates are calculated as second derivative with respect to time;

the one-step prediction of the state may become:

, wherein

，/>，

(5.1.4) propagation ofThe points are as follows:

，

calculating a state predicted value and a one-step prediction error covariance:

ballistic missile feature estimation using 'current' statistical model and methodIs a relationship of (2);

(a) "current" accelerationIn order of right->；

(b) "current" accelerationWhen it is negative, it is added>，

wherein ,is the maximum value of positive acceleration, which is set, +. >Is the set negative acceleration maximum, +.>、/>Acquisition method and->The same is true;

the state predictors are:

，

the one-step prediction error covariance is:

，

wherein Is the white noise sequence covariance;

(5.1.5) measurement update, calculationThe point:

，

，

(5.1.6) propagation volume point:

，

wherein Is an observation function; />Is a measurement state;

(5.1.7) calculating a measurement prediction value:

，

(5.1.8) measurement error covariance matrix and cross covariance matrix:

，

，

wherein Obtaining a sensor error parameter for measuring errors;

(5.1.9) calculating a kalman gain:

，

(5.1.10) calculation of state estimate and estimate covariance:

，

in the above->Is the actual measurement value.

Furthermore, step (5.1) comprises the steps of:

(5.2.1) initializing theParticles initializing the state from a priori distribution>；

(5.2.2) Using a volume Kalman Filter pairThe importance of the individual particles is sampled by the time of the volume Kalman filtering and the iteration of the measuring process, the +.>The time point is marked +.>The particle following proposed distribution function of (a) is: />

, in the formula />Is an importance function; />Mean value of the compliance is +.>Variance is->The two values are respectively the process estimation and covariance estimation of the volume Kalman filtering; / >A measurement that is a volume kalman filter;

(5.2.3) resampling, and obtaining weights of all samples by referring to a particle weight iterative formula:

, in the formula />Is the original particle weight; />Is a state equation;

in order to avoid the problem that probability distribution is too concentrated on the observation sight, adding the speed vector constraint of the moving target, and finally obtaining the double constraint of the observation sight and the speed vector;

is provided withTime particle->Is +.>The particle distance observation line is +.>The particle distance velocity vector line is +.>, and />And->All are gaussian distributions:

, wherein />Extracting and obtaining for observation sight; />；An observation matrix for speed in a state equation; />For solving a point-to-straight line distance function;

the particles areThe weights for the observation line of sight constraints are:

the weights based on the velocity vector constraints are:

due to->And->Two variables that are independent of each other, the total weight based on the observed line-of-sight-velocity vector is therefore:

particles->The updated particle weights are:

, wherein />The original particle weight is not added with a time subscript;

the particles are subjected to weight normalization operation, and then the particles are resampled according to the weight after the normalization step, and the weight normalization operation is performed again Indicating particle->At time->The updated weight of the particle normalization is as follows:

，

(5.2.4) output, calculationThe average number of individual particles based on weight is +.>Information of the time object:

，

(5.2.5) finally obtaining the resamplingThe individual particles are regarded as->The time is input, and an iteration process is achieved through the whole flow of the volume particle filtering; and setting iteration times according to the number of the observation vectors, and obtaining the smooth multistage ballistic missile three-dimensional track points with the same length as the observation vectors after the iteration is finished. />

The beneficial effects are that: compared with the prior art, the invention has the following remarkable advantages: the invention mainly utilizes the observation vector measured by the spaceborne infrared sensor to reconstruct the track of the target according to the straight line reconstruction algorithm, and complements the imperfection of the information; screening suspected tracks by utilizing acceleration characteristic constraint and rationality constraint; compared with the traditional method which needs to collect a large amount of missile prior information to construct a template library, the method effectively reduces the dependence on prior information, and meanwhile, the improved volume particle filtering is used for track smoothing, so that errors can be effectively reduced, and track precision can be improved.

Drawings

FIG. 1 is a schematic flow chart of the present invention.

Description of the embodiments

The technical scheme of the invention is further described below with reference to the accompanying drawings.

As shown in fig. 1, the single-star multistage ballistic missile track reconstruction method based on the straight line reconstruction method comprises the following steps:

the invention discloses a single-star multistage ballistic missile track reconstruction method based on a straight line reconstruction method, which comprises the following steps of:

(1) Converting the target two-dimensional coordinates of the infrared images of the continuous time sequence and the three-dimensional coordinates of the early warning satellite into a unified coordinate system, and extracting a time sequence satellite observation vector;

the step (1) comprises the following steps:

(1.1) converting a target two-dimensional coordinate based on a pixel coordinate system in a satellite image into an image plane coordinate system according to a conversion mode of the pixel coordinate system and the image plane coordinate system, wherein the conversion mode of the two is as follows:

, wherein ,/>Is the position of the target in the image plane coordinate system; />Is the position of the target in the pixel coordinate system; />Is the width of a two-dimensional photo,/>Is the line height of a two-dimensional photograph,/>The size of the pixels of the two-dimensional photo can be obtained from the camera parameters of the satellite;

(1.2) converting the target coordinates in the image plane coordinate system into the satellite body coordinate system according to the linear projection conversion of the image plane coordinate system and the satellite body coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>The position of the target under the satellite body coordinate system; />Representative is the focal length of the sensor, which can be obtained from the camera parameters of the satellite;

(1.3) converting the target coordinates in the satellite body coordinate system into the satellite orbit measurement coordinate system according to the conversion matrix of the satellite body coordinate system and the satellite orbit measurement coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>For the position of the object in the orbit measurement coordinate system,for the position of the object under the body coordinate system, +.>Is based on the attitude angle->The expression mode of the created rotation matrix is as follows: />

, wherein ,/>For observing the pitch angle, the roll angle and the yaw angle of the satellite, the pitch angle, the roll angle and the yaw angle of the satellite can be obtained according to a measurement device built in the satellite;

(1.4) converting the target coordinates in the satellite orbit measurement coordinate system into the WGS84 fixed coordinate system according to a conversion matrix of the satellite orbit measurement coordinate system and the WGS84 fixed coordinate system, wherein the conversion mode is as follows:

, wherein ,/>For satellites in the WGS84 geodetic fixed coordinate system,for the target in the WGS84 fixed coordinate system,/L>Is according to the latitude of the satellite point below>With longitude->The expression mode of the created rotation matrix is as follows:

，

(1.5) set upFor the launch time point of a multi-stage ballistic missile target, < + > >Establishing an observation vector set for a theoretical shutdown point of a target object of the multi-stage ballistic missile: />Wherein is a->Time->Obtained by the formula:

, wherein />Is a proportional coefficient->The orientation element in the camera can be obtained from satellite camera parameters; the observation vector at any time is obtained from the relational expression.

(2) According to the obtained observation vector, estimating the firing angle of the multi-stage ballistic missile according to a firing estimation algorithm;

the step (2) comprises the following steps:

(2.1) coordinates of three-dimensional image points in the WGS84 coordinate systemAnd observation vector->The constructed straight line and radius are +.>Spherical intersection of sphere radius +.>For the earth radius->Estimated altitude with ballistic missile target>And (3) summing; intersection of straight line and sphere->For the three-dimensional coordinate point corresponding to the target emission point, < >>The elements of the formula are as follows:

，

(2.2) calculating according to the formula of the step (2.1) to obtain the intersection pointAbout sphere radius->The functional expression is->Missile targets corresponding to a plurality of observation moments are expressed asObtaining the sphere radius of multiple observation moments according to the formula of the step (2.1)>The functional relation of (2) is +.>；/>

(2.3) parameters due to the orientation of the planeVery small approximately 0, & gt >And (3) withTwo points solve a plane of projection +.>And is directed to plane->The normal vector of (2) is:

；

(2.4) orderAnd->The defined plane of incidence is +.>Then->And->The plane of the shot is defined asThe method comprises the steps of carrying out a first treatment on the surface of the Altering +.in the formula of step (2.3)>2 and->When, i.e. calculate the acquisition plane +.>Normal vector of->And plane->Normal vector of->Plane +.>And plane->Normal vector of->And->Point at the same point and/or>And->Coincide and all pass through the sphere center of the earth and +.>Point, origin of geocentric coordinates +.>The following relation is given:

，

in the formula Is->About->、/>Functional representation of (a) and the sameIs->About->、/>A functional representation of (a);

(2.5) taking three points in the observation process，/>，/>I.e. in the formula of step (2.4)>The value of (2) is 3, and the formula of the step (2.4) is solved by using the Monte Carlo method to obtain +.>，/>，/>Solution of (2); the coordinates of the multi-stage ballistic missile launching point under the WGS84 coordinate system are +.>The normal vector of the coil plane passing through the missile launching point is:

，

in the formula Is equal to +.about.in WGS84 coordinate system>The axes are in the same direction and the direction vector modulo 1, then the shot is:

。

(3) Reconstructing a discrete three-dimensional track of a multi-stage ballistic missile target by adopting a linear reconstruction method according to the satellite observation vector of the continuous time sequence and the calculated shot;

The step (3) comprises the following steps:

(3.1) constructing a camera projection matrix for projecting a multi-stage ballistic missile onto a photographCamera projection matrix->Has the following definition formula:

, wherein />For pixel coordinate system->Scale factor of axis>For pixel coordinate system->Scale factor of axis>As a scale factor, < >>Is the coordinates of the center of the photo in the pixel coordinate system, < >>Coordinates projected in pixel coordinate system for two-dimensional object,/->Is the coordinates of the object in the geodetic fixed coordinate system,

for the coordinate rotation matrix, +.>For coordinate translation matrix, +.>A projection matrix called a camera;

(3.2) assume that along a straight linePoint projection of a mobile multi-stage ballistic missile target +.>And let->Is straight line->In view->Projection onto, due to->Straight line in the image plane +.>On the other hand, the flow of qi is->At this time, a straight line +.>Let view->The expression in the WGS84 solid coordinate system is:

wherein parameter->、/>、/>View +.>Time observation line of sight +.>、/>、；/>Can be substituted into camera coordinates, namely satellite coordinates>Obtaining;

(3.3) setting a straight lineThe parametric equation expression of (2) is:

, wherein />、/>、/>For view->The three-dimensional coordinates of the moving target at the corresponding moment;

along a certain straight line to the targetThree-dimensional position ∈ >And straight line->Direction vector->The following constraint formula is set:

，

，

wherein Normal vector to the warp plane of the target emission point, +.>Is constant (I)>Refers to a priori velocity information of the target; so far, all three-dimensional straight line spaces are embedded in a five-dimensional projection space and are constrained by the constraint formula;

(3.4) passing through a straight lineAnd perpendicular to view->Is>The method comprises the following steps:

thus->The expression of (2) can be solved by:

after rotation of the coordinate system, the straight line is represented under the pixel coordinate system>：

Its parameterization is expressed as +>Is provided with->Is straight line->Go up to any point and letFor->In view->Projection onto; obviously->Because of->Straight line in the image plane +.>Go up and->；

According to camera projection matrix at different momentsAnd a time-series three-dimensional coordinate point of the object->For unknown straight line->The following linear constraint formula is established: />

, wherein />For straight line along the target->Projection on photo, parameterized as +.>The method comprises the steps of carrying out a first treatment on the surface of the When->The unique solution can be solved at this time, and when +.>When the characteristic difference between the actual track and the straight line is too large or more than five views can obtain a least square solution; in general, the least squares estimation matrix has a rank of 5 and +. >To estimate the zero space of the matrix; when->In the time-course of which the first and second contact surfaces,or when the track of the camera projection center is straight during certain tracking periods of the camera movement, the situation that the estimated matrix rank is 4 can occur, and the constraint is acted at the moment; at this time there is a ∈>、/>Extended two-dimensional zero space, thus->Scalar->The optimal value can be calculated according to the constraint formula in the step (3.3);

(3.5) according to straight linesRecursion of parameter equation, reconstructing discrete coordinates of the moving point, and solving straight line according to the linear constraint formula in the step (3.4)>Expressed by the following parametric equation:

，

wherein Is straight line->Go up to a point, let go of>Is straight line->Direction parameters of (2);

straight lineIs>The solution can be found by:

；

(3.6) since the straight line reconstruction algorithm needs to use the information of several frames of continuous time sequence images to accurately reconstruct the straight line track, if the time interval is shorter when the observation satellite starts shooting the ballistic missile target, the straight line reconstruction method is easy to sink into the condition of estimating the matrix lack rank, which leads to more than one result meeting the condition, and all the suspected track results are recorded:

, wherein />Is a suspected trace->Three-dimensional coordinates of the time instant>Is->Time of moment->Taking the value according to the three-dimensional coordinate number >And (5) taking a value according to the number of suspected trajectories.

(4) Constructing acceleration characteristics and rationality limiting conditions of the multi-stage ballistic missile, and screening suspected trajectories according to the acceleration characteristics and rationality limiting conditions of the ballistic missile;

for missiles of the same type, i.e.Different launching programs are adopted, and the ratio of the variable thrust of the missile to the time-varying massAlso has similar characteristic rules, and is usually expressed by acceleration characteristics; therefore, the acceleration characteristics of the multi-stage ballistic missile target are essential characteristics, the ballistic missile also has certain inherent properties and can be constructed into rationality constraint, and suspected trajectory screening is carried out according to the acceleration characteristics and the rationality constraint;

step (4) comprises the following steps:

(4.1) establishing an engine and mass change model, wherein the engine thrust model is as follows:

, wherein ,/>Is the amount of time, representing the moment of movement of the target; />Is the elevation of the target at the current moment; />Is the end time of different active segment stages; />The fuel consumption rate, which is the target of the different active stage phases, and in most cases +.>Is a definite value, can be regarded as->；/>Is the stable exhaust speed of the gas nozzle at different active stage, and can be considered as +. >Is a fixed value; />Is the cross sectional area of the engine nozzle at different active stage; />Is the static pressure of air flow at the nozzles of different active stage; />Is the atmospheric static pressure of the elevation where the target flying process is located; the ballistic missiles of the same model are not very different>、/>、/>、/>、/>The values are only needed to be modified according to a transmitting program in actual early warning, so that the time sequence thrust information of a certain missile can be obtained;

the mass change model is as follows:

, wherein ,/>Is the initial mass of different active section stages, and is generally missile with the same model->、/>The value of (2) is fixed;

(4.2) modeling non-engine acceleration, air resistanceThe general expression of (2) is:

, wherein />Is dynamic pressure (is->Is air density->Is the velocity of the missile relative to air; />Is the maximum cross-sectional area of the missile body; />For Mach number>Record->Is the speed of sound; />As the drag coefficient, the footnote "0" indicates the attack angle of the missile +.>Is 0;

because the included angle between the air resistance and the moving direction of the missile target is a flat angle, the resistance accelerationOpposite to the direction vector of the motion direction, can be expressed as follows:

, wherein ,/>Numerical representation of the characteristics of the missile shape design >As a function of the air density as a function of the elevation, wherein +.>，；/>Is a resistance function, and can represent the force born by the missile according to the speed value;

namely, the analysis and calculation of acceleration caused by the gravity of the earth to the ballistic missile target are realized, the numerical value of the gravity is calculated by using the expression of the space sphere, and the potential energy function obtained by the numerical value is taken asThe method comprises the steps of carrying out a first treatment on the surface of the Let the earth's long half shaft be +.>The short half shaft is +.>And the mass of the earth is symmetrical about the equatorial plane and all longitudinal planes, then the gravitational function of the earth can be expressed as:

, wherein />For the distance between the earth center and the target->For the quality of the object, +.>Gravity inclination angle->For perturbation term is constant, +.>Stress constant->Is Legendre polynomial;

in practice perturbation itemsThe values are not large, so in the earth gravity calculation, get +.>The perturbation term, the gravity function expression at this time, is:

, wherein ,/>；/>

Order theCorresponding gravitational acceleration->Radial in the earth>North gravity dip->Is as follows:

projecting gravitational acceleration to the geocentric sagittal approach>Rotation axis of earthThe following steps are:

；

(4.3) setting active segment ballistic base constraints includes:

outside conditions: the trajectory of a ballistic missile cannot be within the earth, i.e. the elevation of each point on the trajectory needs to be greater than 0, assuming For the elevation of each point on the missile track, the constraint is as follows:

，

elevation increment: the elevation of points on the trajectory of the ballistic missile engine during the working phase must be consistent with the characteristic of increasing with time, and the specific constraint expression is as follows:

，

speed extremum: the maximum scalar rate of points on the missile target trajectory cannot exceed the first cosmic speed because the missile target cannot reach such a high speed with the power of the missile engine, and the first cosmic speed is recorded asThe actual constraint expression is as follows:

，

rate change with time: the scalar rate at each moment on the trajectory of the ballistic missile engine working phase is required to strictly conform to the gradually increasing constraint, and the actual constraint expression is as follows:

acceleration extremum: the maximum scalar acceleration at each moment on the trajectory of the ballistic missile engine working phase cannot be greater than +.>The specific constraint expression is as follows:

，

the initial point elevation condition is the initial detection point elevation condition: due to the characteristics of the infrared sensor, when a cloud layer is shielded, the sensor can hardly detect the tail flame of the ballistic missile; when the sensor can observe the target, the missile generally moves to a high position with thin cloud layer; after azimuth information is obtained according to the observation sight of the initial point, longitude and latitude near the initial point can be calculated, and then elevation of the initial point can be calculated according to the local cloud layer elevation information at the moment Is a constraint range of (2)；

Gas injection limit: the fuel gas injection speed determines the thrust force of the missile engine, and two different forms of fuel injection rate constraints are set, wherein the liquid fuel injection rate is 2500-4600m/s, and the solid injection rate is 2000-3000m/s;

the viewing constraint: the distance between the observation vector and the intersection point of the earth cannot be further than the distance between the missile target and the satellite.

(4.4) the long-distance ballistic missile is the vertical ascending section firstly under most of the flight progress conditions, and the missile can be quickly lifted; then, a program turning section is carried out, and the missile gradually turns to a set target direction at the stage; then the gravity turning section, wherein the thrust vector of the missile is parallel to the movement direction, and the change of the speed vector is only influenced by gravity; finally, the missile is a passive section, and at the moment, the engine of the missile is flamed, and most of external force of the missile comes from the gravity of the earth; in order to better screen suspected trajectories for the multi-stage ballistic missiles, specialized constraints for the ballistic missile targets are set, and the suspected trajectory screening is carried out together with acceleration constraints and rationality constraints according to ballistic missile target rationality constraints consisting of active section ballistic basic constraints and specialized constraints for the ballistic missile targets, so that error trajectory segments can be eliminated;

The specialized constraints for the ballistic missile target in step (4.4) include:

active segment vertical ascent segment durationConstraint:

，

the active section rises vertically Duan GaochengConstraint:

，

first-stage boosting transverse displacementConstraint:

，

duration of the whole active sectionConstraint:

，

whole active section flying elevationConstraint:

，

flight transverse-longitudinal ratio of whole driving sectionConstraint:

。

(5) Performing track smoothing filtering on the discrete multi-stage ballistic missile track by using the volume Kalman filtering with multiple weight judgment to obtain a smooth multi-stage ballistic missile three-dimensional track point with the same length as the observation vector;

step (5) comprises the steps of:

(5.1) transferring the particle sample obtained by prior distribution to a maximum likelihood area by means of an optimal proposal distribution function, thereby completing the aim of optimizing the algorithm performance; the volume particle filtering adopts the state estimation and the estimation covariance of each moment in the volume Kalman filtering step as the mean value and the variance in the particle filtering algorithm to establish a proposal distribution function; the volume Kalman filtering starts to select sigma points according to the estimated state, finally obtains some sigma points with the same weight, and then carries out state estimation and estimation covariance calculation according to the set system model prediction and iteration state value;

The step (5.1) specifically comprises the following steps:

(5.1.1) initializing state and covariance, setting vector length of state equation，/>Linear reconstructed trajectory data for tracking a preliminary stage ballistic missile target>Calculating the state and covariance matrix of the missile:

，/>, wherein />To solve the expectations; />

(5.1.2) time update, calculationPersonal->The point:

，/>, in the formula ：/>Is->Personal->Point (S)>Is thatDecomposing;

(5.1.3) propagationPersonal->Point, one-step prediction of state:

, wherein />In discrete form of a state equation;

inputting a matrix for a process;

and is also provided with，/>Constant determined for the target a priori acceleration means, < +.>Time intervals for satellite observations;

at this time, one-step prediction of the state is regarded as an average value of the state at this time:

, wherein />Expressed as pair->、/>、/>Coordinates are calculated as second derivative with respect to time;

the one-step prediction of the state may become:

, wherein

，/>，

(5.1.4) propagation ofThe points are as follows:

，

calculating a state predicted value and a one-step prediction error covariance:

ballistic missile feature estimation using 'current' statistical model and methodIs a relationship of (2);

(a) "current" accelerationIn order of right->；

(b) "current" accelerationWhen it is negative, it is added>，

wherein ,is the maximum value of positive acceleration, which is set, +. >Is the set negative acceleration maximum, +.>、/>Acquisition method and->The same is true;

the state predictors are:

，

the one-step prediction error covariance is:

，

wherein Is the white noise sequence covariance;

(5.1.5) measurement update, calculationThe point:

，

，

(5.1.6) propagation volume point:

，

wherein Is an observation function; />Is a measurement state;

(5.1.7) calculating a measurement prediction value:

，

(5.1.8) measurement error covariance matrix and cross covariance matrix:

，

，

wherein Obtaining a sensor error parameter for measuring errors;

(5.1.9) calculating a kalman gain:

，

(5.1.10) calculation of state estimate and estimate covariance:

，

in the above->Is the actual measurement value.

(5.2) transferring the particle sample obtained by prior distribution to a maximum likelihood area by means of an optimal proposal distribution function, thereby completing the aim of optimizing the algorithm performance; the volume particle filtering adopts the calculation process of state estimation and estimation covariance in the volume Kalman filtering step (5.1) to obtain state prediction at each moment, and a covariance matrix of the state prediction is used as the mean value and the variance in a particle filtering algorithm to establish a proposal distribution function, and a smooth multistage ballistic missile three-dimensional track point with the same length as an observation vector is obtained after the volume particle filtering iterative calculation;

The step (5.2) specifically comprises the following steps:

(5.2.1) initializing theParticles initializing the state from a priori distribution>；

(5.2.2) Using a volume Kalman Filter pairThe importance of the individual particles is sampled by the time of the volume Kalman filtering and the iteration of the measuring process, the +.>The time point is marked +.>The particle following proposed distribution function of (a) is:

, in the formula />Is an importance function; />Mean value of the compliance is +.>Variance is->The two values are respectively the process estimation and covariance estimation of the volume Kalman filtering; />A measurement that is a volume kalman filter;

(5.2.3) resampling, and obtaining weights of all samples by referring to a particle weight iterative formula:

, in the formula />Is the original particle weight; />Is a state equation;

because of the characteristics of the single star passive observation positioning target, the observation information only has angle information, which means that the observation vector can apply observation sight constraint in the particle weight calculation of the volume Kalman filtering, so that posterior probability distribution represented by particles is more similar to an actual target; in order to avoid the problem that probability distribution is too concentrated on the observation sight, the speed vector constraint of the moving target is added, and finally, the double constraint of the observation sight and the speed vector is obtained;

Is provided withTime particle->Is +.>The particle distance observation line is +.>Granules and particlesThe sub-distance velocity vector line is +.>, and />And->All are gaussian distributions:

, wherein />Extracting and obtaining for observation sight; />；An observation matrix for speed in a state equation; />For solving a point-to-straight line distance function;

the particles areThe weights for the observation line of sight constraints are:

the weights based on the velocity vector constraints are:

due to->And->Are two variables independent of each other, thus baseThe total weight for the observed line-of-sight-velocity vector is:

particles->The updated particle weights are:

, wherein />The original particle weight is not added with a time subscript;

at this time, the weights of the particles are different in size and the absolute difference between the particles is too large, so that the particles need to be subjected to weight normalization operation, and then the particles are resampled according to the weights after the normalization step, and at this time, the particles are still usedIndicating particle->At time->The updated weight of the particle normalization is as follows:

，/>

(5.2.4) output, calculationThe average number of individual particles based on weight is +.>Information of the time object:

，

(5.2.5) finally obtaining the resamplingThe individual particles are regarded as->The time is input, and an iteration process is achieved through the whole flow of the volume particle filtering; and setting iteration times according to the number of the observation vectors, and obtaining the smooth multistage ballistic missile three-dimensional track points with the same length as the observation vectors after the iteration is finished. / >

Claims (10)

1. A single-star multistage trajectory missile track reconstruction method based on a straight line reconstruction method is characterized by comprising the following steps of:

(1) Converting the target two-dimensional coordinates of the infrared images of the continuous time sequence and the three-dimensional coordinates of the early warning satellite into a unified coordinate system, and extracting a time sequence satellite observation vector;

(2) According to the obtained observation vector, estimating the firing angle of the multi-stage ballistic missile according to a firing estimation algorithm;

(3) Reconstructing a discrete three-dimensional track of a multi-stage ballistic missile target by adopting a linear reconstruction method according to the satellite observation vector of the continuous time sequence and the calculated shot;

(4) Constructing acceleration characteristics and rationality limiting conditions of the multi-stage ballistic missile, and screening suspected trajectories according to the acceleration characteristics and rationality limiting conditions of the ballistic missile;

(5) And carrying out track smoothing filtering on the discrete multi-stage ballistic missile track by using the volume Kalman filtering with multiple weight judgment to obtain a smooth multi-stage ballistic missile three-dimensional track point with the same length as the observation vector.

2. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the step (1) comprises the following steps:

(1.1) converting a target two-dimensional coordinate based on a pixel coordinate system in a satellite image into an image plane coordinate system according to a conversion mode of the pixel coordinate system and the image plane coordinate system, wherein the conversion mode of the two is as follows:

, wherein ,/>Is the position of the target in the image plane coordinate system; />Is the position of the target in the pixel coordinate system; />Is the width of a two-dimensional photo,/>Is the line height of a two-dimensional photograph,/>The size of the pixels of the two-dimensional photo can be obtained from the camera parameters of the satellite;

(1.2) converting the target coordinates in the image plane coordinate system into the satellite body coordinate system according to the linear projection conversion of the image plane coordinate system and the satellite body coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>The position of the target under the satellite body coordinate system; />Representative is the focal length of the sensor, which can be taken from the satellite's cameraObtaining the number;

(1.3) converting the target coordinates in the satellite body coordinate system into the satellite orbit measurement coordinate system according to the conversion matrix of the satellite body coordinate system and the satellite orbit measurement coordinate system, wherein the conversion modes of the two are as follows:

, wherein ,/>For the position of the object in the orbit measurement coordinate system, For the position of the object under the body coordinate system, +.>Is based on the attitude angle->The expression mode of the created rotation matrix is as follows:

, wherein ,/>In order to observe the pitch angle, the roll angle and the yaw angle of the satellite, the pitch angle, the roll angle and the yaw angle of the satellite can be obtained according to a measurement device built in the satellite;

(1.4) converting the target coordinates in the satellite orbit measurement coordinate system into the WGS84 fixed coordinate system according to a conversion matrix of the satellite orbit measurement coordinate system and the WGS84 fixed coordinate system, wherein the conversion mode is as follows:

, wherein ,/>For satellites in the WGS84 geodetic fixed coordinate system,for the target in the WGS84 fixed coordinate system,/L>Is according to the latitude of the satellite point below>With longitude->The expression mode of the created rotation matrix is as follows:

(1.5) set->For the launch time point of a multi-stage ballistic missile target, < + >>Establishing an observation vector set for a theoretical shutdown point of a target object of the multi-stage ballistic missile:wherein is a->Time->Obtained by the formula:

, wherein />Is a proportional coefficient->The orientation element in the camera can be obtained from satellite camera parameters; the observation vector at any time is obtained from the relational expression.

3. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the step (2) comprises the following steps:

(2.1) coordinates of three-dimensional image points in the WGS84 coordinate systemAnd observation vector->The constructed straight line and radius are +.>Spherical intersection of sphere radius +.>For the earth radius->Estimated altitude with ballistic missile target>And (3) summing; intersection of straight line and sphere->For the three-dimensional coordinate point corresponding to the target emission point, < >>The elements of the formula are as follows:

(2.2) calculating according to the formula of the step (2.1) to obtain the intersection point +.>About sphere radius->The functional expression is->The missile target corresponding to the observation times is expressed as +.>Obtaining the spherical radius of the multiple observation moments according to the formula of the step (2.1)The functional relation of (2) is +.>；

(2.3) parameters due to the orientation of the planeVery small approximately 0, & gt>And (3) withTwo points solve a plane of projection +.>And is directed to plane->The normal vector of (2) is:

；

(2.4) orderAnd->The defined plane of incidence is +.>Then->And->The defined plane of incidence is +.>The method comprises the steps of carrying out a first treatment on the surface of the Altering +.in the formula of step (2.3)>2 and->When, i.e. calculate the acquisition plane +.>Normal vector of->And plane->Normal vector of (2)Plane +.>And plane->Normal vector of->And->Point at the same point and/or>And->Coincide and all pass through the sphere center of the earth and +.>Point, origin of geocentric coordinates +. >The following relation is given:

，

in the formula Is->About->、/>Functional representation of (a) and the sameIs->About->、/>A functional representation of (a);

(2.5) taking three points in the observation process，/>，/>I.e. in the formula of step (2.4)>The value of (2) is 3, and the formula of the step (2.4) is solved by using the Monte Carlo method to obtain +.>，/>，/>Solution of (2); the coordinates of the multi-stage ballistic missile launching point under the WGS84 coordinate system are +.>The normal vector of the coil plane passing through the missile launching point is:

，

in the formula Is equal to +.about.in WGS84 coordinate system>The axes are in the same direction and the direction vector modulo 1, then the shot is:

。

4. the method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the step (3) comprises the following steps:

(3.1) constructing a camera projection matrix for projecting a multi-stage ballistic missile onto a photographCamera projection matrix->Has the following definition formula:

, wherein />For pixel coordinate system->Scale factor of axis>For pixel coordinate system->Scale factor of axis>As a scale factor, < >>Is the coordinates of the center of the photograph in the pixel coordinate system,coordinates projected in pixel coordinate system for two-dimensional object,/->Is the coordinates of the object in the geodetic fixed coordinate system,

For the coordinate rotation matrix, +.>For coordinate translation matrix, +.>A projection matrix called a camera;

(3.2) assume that along a straight linePoint projection of a mobile multi-stage ballistic missile target +.>And let->Is straight line->In view->Projection onto, due to->Straight line in the image plane +.>On the other hand, the flow of qi is->At this time, a straight line +.>Let viewThe expression in the WGS84 solid coordinate system is:

wherein parameter->、/>、/>View +.>Time observation line of sight +.>、/>、/>；/>Can be substituted into camera coordinates, namely satellite coordinates>Obtaining;

(3.3) setting a straight lineThe parametric equation expression of (2) is:

, wherein />、/>、/>For view->The three-dimensional coordinates of the moving target at the corresponding moment;

along a certain straight line to the targetThree-dimensional position ∈>And straight line->Direction vector->The following constraint formula is set:

，

，

wherein Normal vector to the warp plane of the target emission point, +.>Is constant (I)>Refers to a priori velocity information of the target; so far, all three-dimensional straight line spaces are embedded in a five-dimensional projection space and are constrained by the constraint formula;

(3.4) passing through a straight lineAnd perpendicular to view->Is>The method comprises the following steps:

thus->The expression of (2) can be solved by:

after rotation of the coordinate system, the straight line is represented under the pixel coordinate system >：

Its parameterization is expressed as +>Is provided with->Is straight line->Go up to any point and letFor->In view->Projection onto; obviously->Because of->Straight line in the image plane +.>Go up and->；

According to camera projection matrix at different momentsAnd a time-series three-dimensional coordinate point of the object->For unknown straight line->The following linear constraint formula is established:

, wherein />For straight line along the target->Projection onto a photo, parameterized asThe method comprises the steps of carrying out a first treatment on the surface of the When->The unique solution can be solved at this time, and when +.>When the characteristic difference between the actual track and the straight line is too large or more than five views can obtain a least square solution; in general, the least squares estimation matrix has a rank of 5 and +.>To estimate the zero space of the matrix; when->When the track of the camera projection center is straight during certain tracking periods of the camera motion, the situation that the estimated matrix rank is 4 can occur, and the constraint is acted at the moment; at this time there is a ∈>、/>Extended two-dimensional zero space, thus->Scalar->The optimal value can be obtained by calculation according to the constraint formula in the step (3.3);

(3.5) according to straight linesRecursion of parameter equation, reconstructing discrete coordinates of the moving point, and solving straight line according to the linear constraint formula in the step (3.4) >Expressed by the following parametric equation:

，

wherein Is straight line->Go up to a point, let go of>Is straight line->Direction parameters of (2);

straight lineIs>The solution can be found by:

；

(3.6) recording all suspected track results:

, wherein />Is a suspected trace->The three-dimensional coordinates of the moment in time,is->Time of moment->Taking the value according to the three-dimensional coordinate number>And (5) taking a value according to the number of suspected trajectories.

5. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the step (4) comprises the following steps:

(4.1) establishing an engine and mass change model, wherein the engine thrust model is as follows:

, wherein ,/>Is the amount of time, representing the moment of movement of the target; />Is the elevation of the target at the current moment; />Is the end time of different active segment stages; />The fuel consumption rate, which is the target of the different active stage phases, and in most cases +.>Is a definite value, can be regarded as->；/>Is the stable exhaust speed of the gas nozzle at different active stage, and can be considered as +.>Is a fixed value; />Is the cross sectional area of the engine nozzle at different active stage; />Is the static pressure of air flow at the nozzles of different active stage; / >Is the atmospheric static pressure of the elevation where the target flying process is located; the ballistic missiles of the same model are not very different>、/>、/>、/>、/>The values are only needed to be modified according to a transmitting program in actual early warning, so that the time sequence thrust information of a certain missile can be obtained;

the mass change model is as follows:

, wherein ,/>Is the initial mass of different active section stages, and is generally missile with the same model->、/>The value of (2) is fixed;

(4.2) modeling non-engine acceleration, air resistanceThe general expression of (2) is:

, wherein />Is dynamic pressure (is->Is air density->Is the velocity of the missile relative to air; />Is the maximum cross-sectional area of the missile body; />For Mach number>Record->Is the speed of sound; />As the drag coefficient, the footnote "0" indicates the attack angle of the missile +.>Is 0;

because the included angle between the air resistance and the moving direction of the missile target is a flat angle, the resistance accelerationOpposite to the direction vector of the motion direction, can be expressed as follows:

, wherein ,/>Numerical representation of features originally possessed by a missile shape designThe method comprises the steps of,as a function of the air density as a function of the elevation, wherein +.>，/>；Is a resistance function, and can represent the force born by the missile according to the speed value;

Namely, the analysis and calculation of acceleration caused by the gravity of the earth to the ballistic missile target are realized, the numerical value of the gravity is calculated by using the expression of the space sphere, and the potential energy function obtained by the numerical value is taken asThe method comprises the steps of carrying out a first treatment on the surface of the Let the earth's long half shaft be +.>The short half shaft is +.>And the mass of the earth is symmetrical about the equatorial plane and all longitudinal planes, then the gravitational function of the earth can be expressed as:

, wherein />For the distance between the earth center and the target->For the quality of the object, +.>Gravity inclination angle->For perturbation term is constant, +.>Stress constant->Is Legendre polynomial;

in practice perturbation itemsThe values are not large, so in the earth gravity calculation, get +.>The perturbation term, the gravity function expression at this time, is:

, wherein ,/>；

Order theCorresponding gravitational acceleration->Radial in the earth>North gravity dip->Is as follows:

projecting gravitational acceleration to the geocentric sagittal approach>Earth rotation axis->The following steps are:

；

(4.3) setting an active segment ballistic base constraint,

and (4.4) setting specialized constraint aiming at the ballistic missile target, and carrying out suspected track screening according to the rationality constraint of the ballistic missile target consisting of the initiative section ballistic foundation constraint and the specialized constraint aiming at the ballistic missile target, wherein the acceleration constraint and the rationality constraint together, so that the error track section can be eliminated.

6. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the active segment ballistic basic constraint in step (4.3) comprises:

outside conditions: the trajectory of a ballistic missile cannot be within the earth, i.e. the elevation of each point on the trajectory needs to be greater than 0, assumingFor the elevation of each point on the missile track, the constraint is as follows:

，

elevation increment: the elevation of points on the trajectory of the ballistic missile engine during the working phase must be consistent with the characteristic of increasing with time, and the specific constraint expression is as follows:

，

speed extremum: the maximum scalar rate of points on the missile target trajectory cannot exceed the first cosmic speed because the missile target cannot reach such a high speed with the power of the missile engine, and the first cosmic speed is recorded asThe actual constraint expression is as follows:

，

rate change with time: the scalar rate at each moment on the trajectory of the ballistic missile engine working phase is required to strictly conform to the gradually increasing constraint, and the actual constraint expression is as follows:

acceleration extremum: the maximum scalar acceleration at each moment on the trajectory of the ballistic missile engine working phase cannot be greater than +.>The specific constraint expression is as follows:

，

The initial point elevation condition is the initial detection point elevation condition: after azimuth information is obtained according to the observation sight of the initial point, the longitude and latitude near the initial point can be calculated, and then the elevation of the initial point can be calculated according to the local cloud layer elevation information at the momentConstraint range of->；

Gas injection limit: the fuel gas injection speed determines the thrust force of the missile engine, and two different forms of fuel injection rate constraints are set, wherein the liquid fuel injection rate is 2500-4600m/s, and the solid injection rate is 2000-3000m/s;

the viewing constraint: the distance between the observation vector and the intersection point of the earth cannot be further than the distance between the missile target and the satellite.

7. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the specialized constraints for the ballistic missile target in the step (4.4) include:

active segment vertical ascent segment durationConstraint:

，

the active section rises vertically Duan GaochengConstraint:

，

first-stage boosting transverse displacementConstraint:

，

duration of the whole active sectionConstraint:

，

whole active section flying elevationConstraint:

，

flight transverse-longitudinal ratio of whole driving sectionConstraint:

。

8. the method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: the step (5) comprises the following steps:

(5.1) the volume Kalman filtering starts to select sigma points according to the estimated state, finally obtains some sigma points with the same weight, and then carries out state estimation and estimation covariance calculation according to the set system model prediction and iteration state values;

and (5.2) performing volume particle filtering, namely obtaining state prediction at each moment by adopting a state estimation and covariance estimation calculation process in the volume Kalman filtering step (5.1), taking a covariance matrix of the state prediction as a mean value and a variance in a particle filtering algorithm to establish a proposal distribution function, and obtaining a smooth multistage ballistic missile three-dimensional track point with the same length as an observation vector after volume particle filtering iterative calculation.

9. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: said step (5.1) comprises the steps of:

(5.1.1) initializing state and covariance, setting vector length of state equation，/>Linear reconstructed trajectory data for tracking a preliminary stage ballistic missile target>Calculating the state and covariance matrix of the missile:

，/>, wherein />To solve the expectations;

(5.1.2) time update, calculationPersonal- >The point:

，/>, in the formula ：/>Is->Personal->Point (S)>Is->Decomposing;

(5.1.3) propagationPersonal->Point, one-step prediction of state:

, wherein />In discrete form of a state equation;

inputting a matrix for a process;

and is also provided with，/>A constant value determined for the target a priori acceleration means,time intervals for satellite observations;

at this time, one-step prediction of the state is regarded as an average value of the state at this time:

, wherein />Expressed as pair->、/>、/>Coordinates are calculated as second derivative with respect to time;

the one-step prediction of the state may become:

, wherein

，/>，

(5.1.4) propagation ofThe points are as follows:

，

calculating a state predicted value and a one-step prediction error covariance:

ballistic missile feature estimation using 'current' statistical model and methodIs a relationship of (2);

(a) "current" accelerationIn order of right->；

(b) "current" accelerationWhen it is negative, it is added>，

wherein ,is the maximum value of positive acceleration, which is set, +.>Is the set negative acceleration maximum, +.>、/>Acquisition method and->The same is true;

the state predictors are:

，

the one-step prediction error covariance is:

，

wherein Is the white noise sequence covariance;

(5.1.5) measurement update, calculationThe point:

，

，

(5.1.6) propagation volume point:

，

wherein Is an observation function; />Is a measurement state;

(5.1.7) calculating a measurement prediction value:

，

(5.1.8) measurement error covariance matrix and cross covariance matrix:

，

，

wherein Obtaining a sensor error parameter for measuring errors;

(5.1.9) calculating a kalman gain:

，

(5.1.10) calculation of state estimate and estimate covariance:

，

in the above->Is the actual measurement value.

10. The method for reconstructing the single-star multistage ballistic missile track based on the straight line reconstruction method, which is characterized by comprising the following steps of: said step (5.1) comprises the steps of:

(5.2.1) initializing theParticles initializing the state from a priori distribution>；

(5.2.2) Using a volume Kalman Filter pairThe importance of the individual particles is sampled by the time of the volume Kalman filtering and the iteration of the measuring process, the +.>The time point is marked +.>The particle following proposed distribution function of (a) is:

, in the formula />Is an importance function; />Mean value of the compliance is +.>Variance is->The two values are respectively the process estimation and covariance estimation of the volume Kalman filtering; />A measurement that is a volume kalman filter;

(5.2.3) resampling, and obtaining weights of all samples by referring to a particle weight iterative formula:

, in the formula />Is the original particle weight; />Is a state equation;

in order to avoid the problem that probability distribution is too concentrated on the observation sight, adding the speed vector constraint of the moving target, and finally obtaining the double constraint of the observation sight and the speed vector;

is provided withTime particle->Is +.>The particle distance observation line is +.>The particle distance velocity vector line is +.>, and />And->All are gaussian distributions: />

, wherein />Extracting and obtaining for observation sight; />；/>An observation matrix for speed in a state equation; />For solving a point-to-straight line distance function;

the particles areThe weights for the observation line of sight constraints are:

the weights based on the velocity vector constraints are:

due to->And->Two variables that are independent of each other, the total weight based on the observed line-of-sight-velocity vector is therefore:

particles->The updated particle weights are:

, wherein />The original particle weight is not added with a time subscript;

the particles are subjected to weight normalization operation, and then the particles are resampled according to the weight after the normalization step, and the weight normalization operation is performed againIndicating particle->At time->The updated weight of the particle normalization is as follows:

，

(5.2.4) output, calculationThe average number of individual particles based on weight is +. >Information of the time object:

，

(5.2.5) finally obtaining the resamplingThe individual particles are regarded as->The time is input, and an iteration process is achieved through the whole flow of the volume particle filtering; and setting iteration times according to the number of the observation vectors, and obtaining the smooth multistage ballistic missile three-dimensional track points with the same length as the observation vectors after the iteration is finished. />