CN114326813B - Method and system for predicting residual flight time of unpowered aircraft - Google Patents

Method and system for predicting residual flight time of unpowered aircraft Download PDF

Info

Publication number
CN114326813B
CN114326813B CN202111672205.0A CN202111672205A CN114326813B CN 114326813 B CN114326813 B CN 114326813B CN 202111672205 A CN202111672205 A CN 202111672205A CN 114326813 B CN114326813 B CN 114326813B
Authority
CN
China
Prior art keywords
segment
current segment
flight
calculating
sight
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.)
Active
Application number
CN202111672205.0A
Other languages
Chinese (zh)
Other versions
CN114326813A (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.)
Beijing Aerospace Automatic Control Research Institute
Original Assignee
Beijing Aerospace Automatic Control Research Institute
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 Beijing Aerospace Automatic Control Research Institute filed Critical Beijing Aerospace Automatic Control Research Institute
Priority to CN202111672205.0A priority Critical patent/CN114326813B/en
Publication of CN114326813A publication Critical patent/CN114326813A/en
Application granted granted Critical
Publication of CN114326813B publication Critical patent/CN114326813B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Abstract

The method comprises the steps of obtaining a target distance, a sight inclination angle, a sight deflection angle and a speed vector of the vehicle at the current moment, calculating a total lead angle, and judging turning flight or straight flight of the current flight behavior according to the total lead angle; if the turning flight is performed, segmenting according to the total lead angle, wherein the last segment is a straight line segment, the rest segments are turning segments, and predicting the rest flight time of the turning segments by adopting segmented iteration; calculating the remaining linear segment range, segmenting the linear segment according to the remaining linear segment range, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the turning section residual flight time and the straight line section residual flight time is the predicted total residual flight time; if the straight line flight is adopted, calculating a remaining straight line flight path, segmenting the straight line according to the remaining straight line flight path, and predicting the remaining flight time of the straight line by adopting segmented iteration, wherein the remaining flight time of the straight line is the predicted total remaining flight time.

Description

Method and system for predicting residual flight time of unpowered aircraft
Technical Field
The invention relates to the technical field of flight time control of unpowered aircrafts, in particular to a method and a system for predicting the remaining flight time of an unpowered aircraft.
Background
The aircraft attacks the target at a specified time, which is the basis for implementing a multi-aircraft collaborative saturation attack. To achieve accurate guidance with time constraints, more accurate time of flight remaining information is needed to form the control feedback. Since the remaining time of flight is often difficult to estimate accurately and the speed of an unpowered aircraft varies significantly and uncontrollably in actual flight, in this case, the uncertain change in speed presents certain difficulties in predicting the remaining time of flight. The existing methods predict the remaining flight time based on the assumption that the aircraft speed is constant and the lead angle is small, however, the limitation of the constant speed is too severe and ideal, the flying speed of the unpowered aircraft in the tail section is obviously reduced, and the change process is not controlled. If the method is adopted, the problem of inaccurate estimation of the residual flight time is caused, and the attack time control of the aircraft is failed.
Disclosure of Invention
In view of the above analysis, the embodiments of the present invention aim to provide a method and a system for predicting the remaining time of flight of an unpowered aircraft, which are used for solving the problem of inaccurate estimation of the existing remaining time of flight.
In one aspect, an embodiment of the present invention provides a method for predicting a remaining time of flight of an unpowered aircraft, including the steps of:
acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of an aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight flight of the current flight behavior according to the total lead angle;
if the turning flight is the turning flight, segmenting the turning section according to the total lead angle, and predicting the residual flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
if the straight line flight is performed, calculating a remaining straight line segment course according to the total lead angle and the target distance, segmenting the straight line segment according to the remaining straight line segment course, and predicting the remaining flight time of the straight line segment by adopting segmented iteration, wherein the remaining flight time of the straight line segment is the predicted total remaining flight time.
Based on the further improvement of the technical scheme, the method for predicting the residual flight time of the turning section by adopting the piecewise iteration or predicting the residual flight time of the straight line section by adopting the piecewise iteration comprises the following steps:
for each segment, a state variable of the starting point of the current segment is obtained, wherein the state variable comprises a velocity module value v i Inclination angle theta of sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i Ground level y i
If the current segment is in the turning segment, the state variable further comprises a total lead angle eta i And target distance D tc,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, the state variables further comprise the remaining straight line segment voyage L i
If the current segment is in the turning segment, according to
Figure BDA0003449859010000021
Calculating the initial flight time predicted value of the current segment +.>
Figure BDA0003449859010000022
If the current segment is in the straight line segment, then according to +.>
Figure BDA0003449859010000023
Calculating the initial flight time predicted value of the current segment +.>
Figure BDA0003449859010000024
Wherein Δη i Representing the variation of the total lead angle of the current segment, delta L representing the range variation of the remaining straight line segment of the current segment, K N Representing a navigation ratio;
calculating the initial speed change rate of the current segment according to the dynamics equation of the mass center motion of the aircraft based on the state variable of the starting point of the current segment
Figure BDA0003449859010000025
Based on the initial rate of speed change
Figure BDA0003449859010000031
An initial time-of-flight prediction value for said current segment >
Figure BDA0003449859010000032
Correcting to obtain a current segment flight time correction value +.>
Figure BDA0003449859010000033
State variable based on current segment start point and said current segment time-of-flight correction value
Figure BDA0003449859010000034
Predicting a state variable of a current segment end point, wherein the state variable of the current segment end point is a state variable of a next segment start point;
all segmented time of flight correction values
Figure BDA0003449859010000035
And the remaining time of flight as a predicted turn or straight line segment.
Further, based on the initial rate of speed change
Figure BDA0003449859010000036
An initial time-of-flight prediction value for said current segment>
Figure BDA0003449859010000037
Correcting to obtain a current segment flight time correction value +.>
Figure BDA0003449859010000038
Comprising the following steps:
based on the initial rate of speed change
Figure BDA0003449859010000039
And said initial time of flight prediction value +.>
Figure BDA00034498590100000310
According to the formula
Figure BDA00034498590100000311
Obtaining a first speed predictive value v of the segment end point p, i;
A first speed predictor v based on the segment end point p, i and said initial time of flight prediction value
Figure BDA00034498590100000312
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure BDA00034498590100000313
According to the formula
Figure BDA00034498590100000314
Calculating to obtain a first average correction speed +.>
Figure BDA00034498590100000315
If the current segment is in the turning segment, according to
Figure BDA00034498590100000316
Calculating the current segment time of flight correction +.>
Figure BDA00034498590100000317
If the current segment is in the straight line segment, then according to +.>
Figure BDA00034498590100000318
Calculating the current segment time of flight correction +. >
Figure BDA00034498590100000319
Further, a first speed predictor v based on the segment end point p,i And the initial time-of-flight prediction value
Figure BDA00034498590100000320
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure BDA0003449859010000041
Comprising the following steps:
according to the first speed predictive value v p, i and said initial time of flight prediction value
Figure BDA0003449859010000042
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i
If the current segment is in the turning segment, then according to formula y p,i =-D tc,i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, then according to formula y p,i =-L i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the Wherein D is tc,i+1 Representing the target distance of the starting point of the i+1 segment, L i+1 Representing the rest straight-line segment course of the starting point of the i+1 segment;
based on the first corrected line of sight inclination angle theta Lp,i First corrected line-of-sight offset angle ψ Lp,i And a first speed predictive value v p,i Calculating a first velocity vector in the base reference coordinate systemComponents of each coordinate axis; calculating a first corrected ballistic inclination angle theta based on the components of each coordinate axis of the first velocity vector in the base reference coordinate system p,i
Based on the first speed predictive value v p,i First corrected ground height y p,i First corrected ballistic tilt angle theta p,i Calculating a first corrected speed change rate according to a kinetic equation of the mass center motion of the aircraft
Figure BDA0003449859010000043
Further, according to the first speed predictive value v p,i And the initial time-of-flight prediction value
Figure BDA0003449859010000044
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i Comprising:
correcting the sight inclination angle of the current segment according to the following formula to obtain a first corrected sight inclination angle theta Lp,i
If the current segment is in the turning segment, according to
Figure BDA0003449859010000045
Calculating the initial gaze tilt rate +.>
Figure BDA0003449859010000051
And a first corrected gaze inclination rate of change +>
Figure BDA0003449859010000052
If the current segment is in the straight line segment, then according to
Figure BDA0003449859010000053
Calculating the initial gaze tilt rate +.>
Figure BDA0003449859010000054
And a first corrected gaze inclination rate of change +>
Figure BDA0003449859010000055
According to the formula
Figure BDA0003449859010000056
Obtaining a first corrected sight inclination angle theta Lp,i
Correcting the sight line deflection angle of the current segment according to the following formula to obtain a first corrected sight line deflection angle phi Lp,i
If the current segment is in the turning segment, according to
Figure BDA0003449859010000057
Figure BDA0003449859010000058
Calculating the initial line of sight deflection change rate +.>
Figure BDA0003449859010000059
And a first corrected line-of-sight angle change rate
Figure BDA00034498590100000510
If the current segment is in the straight line segment, then according to
Figure BDA00034498590100000511
Figure BDA00034498590100000512
Calculating the initial line of sight deflection change rate +.>
Figure BDA00034498590100000513
And a first corrected line-of-sight angle change rate
Figure BDA00034498590100000514
According to the formula
Figure BDA00034498590100000515
Obtaining a first corrected sight offset angle psi Lp,i
Further, a state variable based on a current segment start point and the current segment time of flight correction value
Figure BDA00034498590100000516
Predicting a state variable for a current segment end point, comprising:
initial rate of speed change according to current segment
Figure BDA0003449859010000061
And said current segment time of flight correction value +.>
Figure BDA0003449859010000062
According to the formula->
Figure BDA0003449859010000063
Obtaining a second speed predicted value v of the current segment end point q,i
A second velocity predictor v based on the segment end point q,i And the current segment time of flight correction value
Figure BDA0003449859010000064
Calculating a second rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure BDA0003449859010000065
For initial rate of speed change
Figure BDA0003449859010000066
And a second corrected speed change rate->
Figure BDA0003449859010000067
Averaging to obtain a second average correction speed
Figure BDA0003449859010000068
Based on the second average corrected speed
Figure BDA0003449859010000069
Calculating a third speed predicted value of the end point of the current segment
Figure BDA00034498590100000610
The third speed predicted value is the predicted speed v of the ending point of the current segment i+1
Line of sight inclination angle theta according to current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value
Figure BDA00034498590100000611
Calculating the sight inclination angle theta of the end point of the current segment L,i+1 The method comprises the steps of carrying out a first treatment on the surface of the Line-of-sight offset angle psi according to the current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value +.>
Figure BDA00034498590100000612
Calculating to obtain the sight offset angle psi of the end point of the current segment L,i +1;
If the current segment is in the turning segment, then according to formula y i+1 =-D tc,i+1 sinθ L,i+1 Calculating the ground height of the ending point of the current segment; if the current segment is in a straight line segment, then according to y i+1 =-L i+1 sinθ L,i+1 The ground height of the end point of the current segment is calculated.
Further, based on the state variable of the starting point of the current segment, calculating the initial speed change rate of the current segment according to the dynamics equation of the movement of the mass center of the aircraft
Figure BDA00034498590100000613
Comprising the following steps:
based on the inclination angle theta of the sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i And velocity module v i Calculating components of the speed vector in each coordinate axis in the base reference coordinate system; calculating the component of each coordinate axis in the basic reference coordinate system based on the velocity vector to obtain the initial ballistic inclination angle theta of the current segment i
If the current segment is in the turning segment, the method is based on the formula
Figure BDA0003449859010000071
Figure BDA0003449859010000072
Calculating components of acceleration in a ballistic coordinate system; if the current segment is in the straight line segment, then according to the formula +.>
Figure BDA0003449859010000073
Calculating components of acceleration in a ballistic coordinate system;
according to the ground height y i Interpolation to obtain the atmospheric density ρ and the sound speed aT, and Mach number is calculated according to the atmospheric density ρ and the sound speed aT
Figure BDA0003449859010000074
And dynamic pressure->
Figure BDA0003449859010000075
Calculating the lift coefficient C according to the following formula Ld,i
F y =MT·g tc cosθ i +MT·a y,i ,F z =MT·a z,i
Figure BDA0003449859010000076
Figure BDA0003449859010000077
Based on Mach and lift coefficient C Ld,i Calculating the attack angle command alpha by Newton iteration method c,i
Based on the Mach number Mach and the angle of attack command alpha c,i Calculating a resistance coefficient C according to the fitting formula D,i
According to the formula
Figure BDA0003449859010000078
The rate of change of the initial velocity is calculated,
wherein MT represents the mass of the aircraft, ST represents the reference area of the aircraft, g tc Representing the gravitational constant.
Further, after the state variable of the current segment start point is obtained for each segment, the initial speed change rate is based on
Figure BDA0003449859010000079
An initial time-of-flight prediction value for said current segment>
Figure BDA00034498590100000710
Before the correction, the method further comprises the following steps:
if the current segment is in the turning segment, calculating the prepositioning angle eta of the ending point of the current segment according to the following formula z,i+1 Front deflection angle eta y,i+1 And a total lead angle eta i+1
Figure BDA0003449859010000081
Figure BDA0003449859010000082
η y,i+1 =η y,i +Δη y,i ,η z,i+1 =η z,i +Δη z,i ,η i+1 =arccos(cosη y,i+1 cosη z,i+1 );
According to the formula
Figure BDA0003449859010000083
Calculating a target distance D of the end point of the current segment tc,i+1
If the current segment is in the straight line segment, calculating the remaining straight line segment range L of the ending point of the current segment according to the following formula i+1 Front tilt angle eta z,i+1 And a prepositive deflection angle eta y,i+1
L i+1 =L i -ΔL,
Figure BDA0003449859010000084
Wherein K is N Represents the navigation ratio, delta η Indicating the pre-angle segment quantity constant and NumP indicates the number of segments of the turn segment.
Further, calculating a total lead angle based on the gaze inclination angle, gaze deflection angle, and velocity vector, comprising:
Calculating a direction cosine matrix from the reference coordinate system to the sight line coordinate system according to the sight line inclination angle and the sight line deflection angle;
calculating a component of the velocity vector in the line-of-sight coordinate system based on the direction cosine matrix;
the total lead angle is calculated according to the following formula:
Figure BDA0003449859010000085
η=arccos(cosη y cosη z )
wherein v is xL 、v yL And v zL Respectively representing the components of the velocity vector in three coordinate axes of a sight line coordinate system, eta z Represents the prepositive dip angle eta y Represents the lead angle, and η represents the total lead angle.
Compared with the prior art, the method and the device have the advantages that the current flying behavior is judged to turn or fly straight according to the total lead angle according to the three-dimensional mass center equation of motion of the aircraft in a closed loop mode, the current flying behavior is segmented according to the total lead angle when the aircraft turns, the remaining flying time is segmented according to the remaining straight-line segment voyage when the aircraft turns, the remaining flying time is predicted by adopting segmented iteration, the future speed of the unpowered aircraft is predicted in a segmented mode through iteration solution, and the remaining flying time estimation formula is corrected in an iteration mode, so that the estimation of the remaining flying time is more accurate, and the remaining flying time can still be accurately predicted under the conditions that the total lead angle greatly changes and the flying speed of the unpowered aircraft is not controlled.
In another aspect, an embodiment of the present invention provides a system for predicting a remaining time of flight of an unpowered aircraft, including:
the data acquisition module is used for acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of the aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight line flight of the current flight behavior according to the total lead angle;
the turning flight prediction module is used for segmenting the turning section according to the total front angle if turning flight is performed, and predicting the remaining flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
and the linear flight prediction module is used for calculating a remaining linear segment course according to the total lead angle and the target distance if the linear flight is linear flight, segmenting the linear segment according to the remaining linear segment course, and predicting the remaining flight time of the linear segment by adopting segmentation iteration, wherein the remaining flight time of the linear segment is the predicted total remaining flight time.
In the invention, the technical schemes can be mutually combined to realize more preferable combination schemes. Additional features and advantages of the invention will be set forth in the description which follows, and in part will be obvious from the description, or may be learned by practice of the invention. The objectives and other advantages of the invention may be realized and attained by the structure particularly pointed out in the written description and drawings.
Drawings
The drawings are only for purposes of illustrating particular embodiments and are not to be construed as limiting the invention, like reference numerals being used to refer to like parts throughout the several views.
FIG. 1 is a flow chart of a method of predicting remaining time of flight of an unpowered aircraft in accordance with an embodiment of the present invention;
FIG. 2 is a block diagram of a residual time of flight prediction system for an unpowered aircraft in accordance with an embodiment of the present invention;
FIG. 3 is a coordinate system diagram of an embodiment of the present invention;
fig. 4 is a graph showing comparison of results of different residual time estimation methods.
Detailed Description
Preferred embodiments of the present invention will now be described in detail with reference to the accompanying drawings, which form a part hereof, and together with the description serve to explain the principles of the invention, and are not intended to limit the scope of the invention.
In order to solve the problem that the prediction of the residual time of the unpowered aircraft is inaccurate under the conditions that the speed changes along with the flight state and the leading angle is large, the current flying behavior turning flight or the straight line flight is judged according to the total leading angle according to a three-dimensional mass center equation of motion of the aircraft in a closed loop mode, the current flying behavior turning flight or the straight line flight is segmented according to the total leading angle when the current flying behavior turning flight is turning flight, the residual flight time is segmented according to the residual straight line section range when the current flying behavior turning flight is straight line flight, the segmentation iteration prediction is adopted, the segmentation prediction is carried out on the future speed of the unpowered aircraft through iteration solution, and the residual flight time estimation formula is iteratively corrected, so that the estimation of the residual flight time is more accurate, and the residual flight time can still be accurately predicted under the conditions that the total leading angle greatly changes and the flight speed of the unpowered aircraft is uncontrolled.
As shown in fig. 3, the coordinate system of the embodiment of the present invention includes a base reference coordinate system, a line-of-sight coordinate system, and a trajectory coordinate system.
The three-dimensional centroid of the aircraft is taken as an origin O in a base reference coordinate system, the X-axis direction is parallel to the horizontal plane and the pointing north is positive, and the base reference coordinate system is taken as X R The Y axis is vertically upward, in Y R The Z axis is parallel to the horizontal plane and is oriented according to the right hand rule, and Z is used for R And (3) representing.
The line-of-sight coordinate system takes the three-dimensional centroid of the aircraft as an origin O, the X-axis direction points to the target point from the aircraft, and the X L Indicating Z L The axial direction being parallel to the horizontal plane and pointing at X L Right is positive, the Y axis determines the direction according to the right hand rule, and Y is used for L And (3) representing.
The trajectory coordinate system takes the three-dimensional centroid of the aircraft as an origin O, and the X-axis direction is along the direction of a velocity vector, and takes X as the direction of the velocity vector V Representing Y V The axis contains the velocity vector and is vertical to the velocity vector and positive upwards in the vertical plane, the Z axis determines the direction according to the right hand rule, and Z is used for V And (3) representing.
Symbols used in the embodiments of the present application are collectively described as follows: d (D) tc Representing a linear distance between the aircraft and the target; f represents stress; f (F) L Representing the total lift; c (C) L Representing the lift coefficient; c (C) D Representing the resistance coefficient; η represents the total lead angle, i.e. the aircraft velocity vector and X L Is included in the plane of the first part; η (eta) y Representing the angular offset (velocity vector in line of sight coordinate system X L -O-Z L In-plane projection and X L An angle between the axes); η (eta) z Representing the pretilt angle (velocity vector and line-of-sight coordinate system X L -O-Z L Included angle of plane); Δη represents the total lead angle variation; Δη y Representing the front deflection angle variation; Δη z Representing the change amount of the pre-tilt angle; delta η A segment quantity constant for the lead angle; v represents a velocity module;
Figure BDA0003449859010000111
Representing the rate of change of speed; />
Figure BDA0003449859010000112
Representing the average rate of change of speed; y represents the height from the ground; θ L Representing inclination of line of sight, i.e. line of sight vector X L And X is R -O-Z R An included angle of the planes; psi phi type L Representing the angle of gaze deviation, i.e. gaze vector X L At X R -O-Z R In-plane projection and X R An included angle of the shaft; />
Figure BDA0003449859010000113
And->
Figure BDA0003449859010000114
The change rate of the inclination angle of the line of sight and the change rate of the deflection angle of the line of sight are respectively represented; θ represents the ballistic dip, i.e., the angle of the velocity vector with the horizontal; a, a y And a z Respectively represent the acceleration Y in a ballistic coordinate system v Axis and Z v The component of the axis, L represents the remaining linear segment range, ΔL represents the linear segment segmentation constant; the index i indicates the current segment and index i+1 indicates the next segment.
In one embodiment of the present invention, a method for predicting the remaining time of flight of an unpowered aircraft is disclosed, as shown in fig. 1, comprising the steps of:
s1, acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of an aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight flight of the current flight behavior according to the total lead angle.
S2, if turning flight is performed, segmenting the turning section according to the total lead angle, and predicting the residual flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
S3, if the straight line flight is performed, calculating a remaining straight line section course according to the total lead angle and the target distance, segmenting the straight line section according to the remaining straight line section course, and predicting the remaining flight time of the straight line section by adopting segmented iteration, wherein the remaining flight time of the straight line section is the predicted total remaining flight time.
The current flying behavior is judged to turn or fly straight according to the total lead angle, the current flying behavior is segmented according to the total lead angle when turning, the current flying behavior is segmented according to the residual straight line range when flying straight, the increment of the lead inclination angle of the aircraft in each segmented interval is guaranteed to be a small angle, and the residual flight time is predicted by adopting segmentation iteration, so that the residual flight time can still be accurately predicted under the conditions that the total lead angle is greatly changed and the flying speed of the unpowered aircraft is not controlled.
Wherein the acquired velocity vector is a velocity vector in a reference base coordinate system in which the components of the three coordinate axes can be represented as v, respectively xR 、v yR And v zR
Specifically, in step S1, calculating the total lead angle based on the gaze inclination angle, the gaze offset angle, and the velocity vector includes:
s11, calculating a direction cosine matrix from the reference coordinate system to the sight line coordinate system according to the sight line inclination angle and the sight line deflection angle;
Specifically, according to the formula
Figure BDA0003449859010000131
Direction cosine matrix for calculating base reference coordinate system to sight line coordinate system>
Figure BDA0003449859010000135
S12, calculating the components of the speed vector in the sight line coordinate system based on the direction cosine matrix;
the components of the velocity vector in the line-of-sight coordinate system are
Figure BDA0003449859010000132
Wherein v is xL 、v yL And v zL Representing the coordinate components of the velocity vector in the line of sight coordinate system, respectively.
S13, calculating a total lead angle according to the following formula:
Figure BDA0003449859010000133
η=arccos(cosη y cosη z )
wherein v is xL 、v yL And v zL Respectively representing the components of the velocity vector in three coordinate axes of a sight line coordinate system, eta z Represents the prepositive dip angle eta y Represents the lead angle, and η represents the total lead angle.
In practice, the pre-angle segment quantity constant may take a smaller quantity, e.g., Δ η 2 ° may be taken. When the total lead angle eta 0 >Δ η When the vehicle is in the turning flight, judging that the current flight is the turning flight, namely, the vehicle is in the turning flight at the current predicted moment of the residual flight; otherwise, judging the current flying behavior to fly linearly.
If the current flight is turning flight, segmenting the turning section according to the total lead angle according to the formula
Figure BDA0003449859010000134
The number of segments of the turning section is calculated as NumP, and floor (·) represents a rounding down. And predicting the flight time of each turn section in sequence to obtain the residual flight time of the turn section, then calculating the range of the residual straight line section according to the total lead angle and the target distance of the starting point of the straight line section, segmenting the straight line section according to the range of the residual straight line section, and predicting the flight time of each straight line section in sequence to obtain the residual flight time of the straight line section, wherein the total residual flight time is the sum of the residual flight time of the turn section and the residual flight time of the straight line section. Wherein the starting point of the straight line segment is the ending point of the turning segment.
If the current flight is straight-line flight, the straight-line segment is segmented directly according to the remaining straight-line segment range. According to the formula
Figure BDA0003449859010000141
Calculating the remaining straight line range according to the formula +.>
Figure BDA0003449859010000142
Calculating the segmentation number NumL+1 of the straight line segment, and predicting the flight time for each segment in sequence to obtain the residual flight time of the straight line segment, namely the predicted total residual flight time.
Wherein D is tc Aircraft representing predicted moments of remaining flightThe linear distance between the targets, η, represents the total lead angle at which the remaining flight prediction is performed, K N For the navigational ratio of the proportional navigational method employed, L represents the remaining straight-line segment voyage at the time of the remaining flight prediction.
Specifically, predicting the residual flight time of the turning section by adopting piecewise iteration or predicting the residual flight time of the straight line section by adopting piecewise iteration comprises the following steps:
s21, for each segment, acquiring a state variable of the starting point of the current segment, wherein the state variable comprises a speed module value v i Inclination angle theta of sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i Ground level y i
If the current segment is in the turning segment, the state variable further comprises a total lead angle eta i And target distance D tc,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, the state variables further comprise the remaining straight line segment voyage L i
If the current flight is turning flight, for the first segment of the turning segment, the state variable of the starting point of the first segment of the straight line segment is the state variable of the ending point of the last segment of the turning segment, wherein the state variable of the starting point of the first segment of the straight line segment is the state value of the predicted moment of the remaining flight
Figure BDA0003449859010000143
D tc,Nump+1 And eta Nump+1 Is the distance and total lead angle of the target at the end point of the last leg of the turn.
And if the current flight is straight-line flight, the state variable of the first segmentation starting point of the straight-line segment is the state value of the rest flight prediction moment.
S22, if the current segment is in the turning segment, according to
Figure BDA0003449859010000151
Calculating the initial flight time predicted value of the current segment +.>
Figure BDA0003449859010000152
If the current segment is in the straight line segment, then according to +.>
Figure BDA0003449859010000153
Calculating the initial flight time predicted value of the current segment +.>
Figure BDA0003449859010000154
Wherein Δη i Representing the variation of the total lead angle of the current segment, delta L representing the range variation of the remaining straight line segment of the current segment, K N Representing a navigation ratio;
specifically, if the current segment is in the turning segment, the change amount Δη of the total lead angle of the current segment i Calculated according to the following formula
Figure BDA0003449859010000155
I.e. when i=nump, i.e. NumP segment takes Δη i =-(η i -NumP·Δ η ) In the 1 st to NumP-1 th paragraphs, Δη is taken i =-Δ η
Figure BDA0003449859010000156
The method is obtained according to the following steps:
para eta i =arccos(cosη y,i cosη z,i ) Differentiation is performed to obtain
Figure BDA0003449859010000157
And because of
Figure BDA0003449859010000158
Figure BDA0003449859010000159
Figure BDA00034498590100001510
Has the following components
Figure BDA0003449859010000161
The movement of the aircraft in three-dimensional space thus has a movement pattern of
Figure BDA0003449859010000162
The elimination time variable is +.>
Figure BDA0003449859010000163
Then in interval t i ,t]Upper integration and reduction->
Figure BDA0003449859010000164
Therefore there is->
Figure BDA0003449859010000165
Then to sin eta t One-order Taylor expansion, sin η t =sinη i +Δη t ·cosη t Substituted with +.>
Figure BDA0003449859010000166
Then at->
Figure BDA0003449859010000167
Upper integral is available->
Figure BDA0003449859010000168
According to the formula
Figure BDA0003449859010000169
Predicting the flight time of the current segment to obtain an initial flight time predicted value +.>
Figure BDA00034498590100001610
If the current segment is in the straight line segment, then according to
Figure BDA00034498590100001611
Predicting the time of flight of the pre-segment to obtain an initial time of flight prediction value for the current segment>
Figure BDA00034498590100001612
S23, calculating the initial speed change rate of the current segment according to the dynamics equation of the movement of the mass center of the aircraft based on the state variable of the starting point of the current segment
Figure BDA00034498590100001613
S24, based on the initial speed change rate
Figure BDA00034498590100001614
An initial time-of-flight prediction value for said current segment>
Figure BDA0003449859010000171
Correcting to obtain a current segment flight time correction value +.>
Figure BDA0003449859010000172
S25, based on the state variable of the current segment start point and the current segment flight time correction value
Figure BDA0003449859010000173
Predicting a state variable of a current segment end point, wherein the state variable of the current segment end point is a state variable of a next segment start point;
If the current segment is the last segment of the straight line segment, the corrected value of the flight time of the current segment is obtained
Figure BDA0003449859010000174
The total remaining time of flight is obtained and the state variable of the current segment end point is no longer predicted. />
All segmented time of flight correction values
Figure BDA0003449859010000175
And the remaining time of flight as a predicted turn or straight line segment.
After step S21, before step S24, further comprising:
s201, if the current segment is in the turning segment, calculating the prepositioning inclination angle eta of the end point of the current segment according to the following formula z,i+1 Front deflection angle eta y,i+1 And a total lead angle eta i+1
Figure BDA0003449859010000176
Figure BDA0003449859010000177
η y,i+1 =η y,i +Δη y,i ,η z,i+1 =η z,i +Δη z,i ,η i+1 =arccos(cosη y,i+1 cosη z,i+1 );
According to the formula
Figure BDA0003449859010000178
Calculating a target distance D of the end point of the current segment tc,i+1
If the current segment is in the straight line segment, calculating the remaining straight line segment range L of the ending point of the current segment according to the following formula i+1 Front tilt angle eta z,i+1 And a prepositive deflection angle eta y,i+1
L i+1 =L i -ΔL,
Figure BDA0003449859010000179
Wherein K is N Represents the navigation ratio, delta η Representing the lead angle segment quantity constant.
In step S23, based on the state variable of the starting point of the current segment, the initial speed change rate of the current segment is calculated according to the dynamics equation of the centroid motion of the aircraft, including steps S231-S237:
s231 based on the sight inclination angle theta L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i And velocity module v i Calculating components of the speed vector in each coordinate axis in the base reference coordinate system; calculating the component of each coordinate axis in the basic reference coordinate system based on the velocity vector to obtain the ballistic inclination angle theta of the current segment i
Specifically, the ballistic tilt θ is calculated according to the following formula:
Figure BDA0003449859010000181
Figure BDA0003449859010000182
wherein, the liquid crystal display device comprises a liquid crystal display device,
Figure BDA0003449859010000183
direction cosine matrix representing base reference coordinate system to line-of-sight coordinate system, v represents velocity module value, eta z Represents the prepositive dip angle eta z Represents the prepositive deflection angle theta L Represents the inclination of the sight line, ψ L Representing the angle of deflection of the line of sight, v px ,v py ,v pz Representing the component sum v of the velocity vector on each axis of the line-of-sight coordinate system qx ,v qy ,v qz Representing the components of the velocity vector on each axis of the reference frame, respectively.
The above formula for calculating the ballistic tilt angle θ is expressed as θ=f θLLyz ,v)。
The state variable line-of-sight inclination angle theta according to the current segment start point L,i Sight line deflection angle psi L,i Velocity module v i Front tilt angle eta z,i And a prepositive deflection angle eta y,i Carry-in θ=f θLLyz V) the initial trajectory of the current segment can be calculatedInclination angle theta i
S232, if the current segment is in the turning segment, according to the formula
Figure BDA0003449859010000184
Figure BDA0003449859010000191
Calculating the component of the current acceleration in a ballistic coordinate system; if the current segment is in the straight line segment, then according to the formula +.>
Figure BDA0003449859010000192
Calculating the component of the current acceleration in a ballistic coordinate system;
a y Y representing acceleration in ballistic coordinate system v Component of axis, a z Z representing acceleration in ballistic coordinate system v Component of the axis.
S233, according to the ground height y i Interpolation to obtain the atmospheric density ρ and the sound speed aT, and Mach number is calculated according to the atmospheric density ρ and the sound speed aT
Figure BDA0003449859010000193
And dynamic pressure->
Figure BDA0003449859010000194
Specifically, according to the U.S. standard atmosphere model, the atmosphere density ρ and the sound speed aT can be obtained by interpolating according to the height of the current segment from the ground
Figure BDA0003449859010000195
Mach number is calculated according to +.>
Figure BDA0003449859010000196
Dynamic pressure is calculated.
S234, calculating a lift coefficient C according to the following formula Ld,i
F y =MT·g tc cosθ i +MT·a y,i ,F z =MT·a z,i
Figure BDA0003449859010000197
Wherein MT represents aircraft mass, ST represents aircraft reference area, g tc Represents the gravitational constant, sign represents the sign function, F L Representing total lift force, F y Y representing lift in ballistic coordinate system v Component of axis, F z Z representing lift in ballistic coordinate system v A component of the shaft; according to the calculation principle of lift force: lift = lift coefficient x dynamic pressure x aircraft reference area, calculate lift coefficient C Ld,i
S235, based on Mach number and lift coefficient C Ld,i Calculating the attack angle command alpha by Newton iteration method c,i
Specifically, the iteration formula of the newton iteration method is as follows:
Figure BDA0003449859010000201
the iteration initial value is 5 DEG, namely alpha 0 At 5 deg., the iteration accuracy meets the accuracy requirement, e.g. 1e-5, or a maximum number of iterations is reached, e.g. 100, the iteration is ended.
Specifically, the function f (α) and derivative f '(α) used by newton's iterative method are as follows:
f(α)=k Cl_00 +k Cl_10 α+k Cl_01 Mach+k Cl_20 α 2 +k Cl_11 α·Mach+k Cl_02 Mach 2 +k Cl_30 α 3
+k Cl_21 α 2 ·Mach+k Cl_12 α·Mach 2 +k Cl_03 Mach 3 +k Cl_31 α 3 ·Mach+k Cl_22 α 2 ·Mach 2
+k Cl_13 α·Mach 3 +k Cl_04 Mach 4 +k Cl_32 α 3 ·Mach 2 +k Cl_23 α 2 ·Mach 3
+k Cl_14 α·Mach 4 +k Cl_05 Mach 5 -C Ld
f'(α)=k Cf_00 +k Cf_10 α+k Cf_01 Mach+k Cf_20 α 2 +k Cf_11 α·Mach+k Cf_02 Mach 2
+k Cf_21 α 2 ·Mach+k Cf_12 α·Mach 2 +k Cf_03 Mach 3 +k Cf_22 α 2 ·Mach 2
+k Cf_13 α·Mach 3 +k Cf_04 Mach 4
wherein k is Cl_00 、k Cl_10 、k Cl_01 、k Cl_20 、k Cl_11 、k Cl_02 、k Cl_30 、k Cl_21 、k Cl_12 、k Cl_03 、k Cl_31 、k Cl_22 、k Cl_13 、k Cl_04 、k Cl_32 、k Cl_23 、k Cl_14 、k Cl_05 ,k Cf_00 、k Cf_10 、k Cf_01 、k Cf_20 、k Cf_11 、k Cf_02 、k Cf_21 、k Cf_12 、k Cf_03 、k Cf_22 、k Cf_13 、k Cf_04 Is a coefficient fitted in advance according to aerodynamic parameters of the aircraft.
Mach number Mach calculated in step S232 and lift coefficient C calculated in step S234 Ld,i Carrying out iterative calculation to obtain attack angle command alpha c,i . The Newton iteration method is high in convergence speed, and a more accurate numerical solution can be obtained by setting convergence conditions, so that a common attack angle instruction can be calculated rapidly and accurately.
S236, based on Mach number and angle of attack instruction alpha c,i Calculating a resistance coefficient C according to the fitting formula D,i
Specifically, the drag coefficient is calculated according to the following formula:
Figure BDA0003449859010000202
wherein k is tc_00 、k tc_10 、k tc_01 、k tc_20 、k tc_11 、k tc_02 、k tc_30 、k tc_21 、k tc_12 、k tc_03 、k tc_40 、k tc_31 、k tc_22 、k tc_13 、k tc_04 、k tc_50 、k tc_ 41、k tc_32 、k tc_23 、k tc_14 、k tc_05 Is a coefficient fitted in advance according to aerodynamic parameters of the aircraft. The aerodynamic characteristics of the aircraft are generally described by a large amount of discrete data, and by fitting the data, an approximate fitting relation between the drag coefficient and the attack angle and Mach number can be obtained, so that the method is simple, and the calculation speed is high.
Mach number Mach calculated in step S232 and attack angle command alpha calculated in step S236 c, i carrying out iterative calculation to obtain a resistance coefficient C D,i
S237, according to the formula
Figure BDA0003449859010000211
The initial rate of speed change is calculated.
In particular, according to the kinetic equation of the movement of the aircraft mass centre
Figure BDA0003449859010000212
And calculating to obtain the initial speed change rate. The stress of the unpowered aircraft mainly comes from resistance and gravity, so the dynamic equation is expressed as
Figure BDA0003449859010000213
Wherein->
Figure BDA0003449859010000214
Represents the effect of resistance on the rate of change of speed g tc sinθ i The effect of gravity on the rate of change of speed is shown.
Wherein MT represents the mass of the aircraft, ST represents the reference area of the aircraft, g tc Representing the gravitational constant.
The method of calculating the rate of change of speed in steps 233-237 is represented as
Figure BDA0003449859010000215
v denotes the velocity modulus in the base reference frame, a y Y representing acceleration in ballistic coordinate system v Component of axis, a z Z representing acceleration in ballistic coordinate system v The component of the axis, θ, represents the ballistic dip and y represents the ground level.
Specifically, step S24 is based on the initial speed change rate
Figure BDA0003449859010000216
Initial time-of-flight predictions for the current segment
Figure BDA0003449859010000217
Correcting to obtain a current segment flight time correction value +.>
Figure BDA0003449859010000218
Including steps S241-S244.
S241, based on the initial speed change rate
Figure BDA0003449859010000221
And said initial time of flight prediction value +.>
Figure BDA0003449859010000222
According to the formula
Figure BDA0003449859010000223
Obtaining a first speed predictive value v of the segment end point p,i
S242, predicting value v of first speed based on the segment ending point p,i And the initial time-of-flight prediction value
Figure BDA0003449859010000224
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid >
Figure BDA0003449859010000225
/>
Specifically, in step S242, a first velocity prediction value v based on the segment end point p,i And the initial time-of-flight prediction value
Figure BDA0003449859010000226
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure BDA0003449859010000227
Comprises steps S2421-S2424:
s2421 based on the first speed predictor v p,i And the initial time-of-flight prediction value
Figure BDA0003449859010000228
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i
Specifically, according to the first speed predictor v p,i And the initial time-of-flight prediction value
Figure BDA0003449859010000229
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i Comprising:
correcting the sight inclination angle of the current segment according to the following formula to obtain a first corrected sight inclination angle theta Lp,i
If the current segment is in the turning segment, according to
Figure BDA00034498590100002210
Calculating the initial gaze tilt rate +.>
Figure BDA00034498590100002211
And a first corrected line of sightInclination angle change Rate->
Figure BDA00034498590100002212
If the current segment is in the straight line segment, then according to
Figure BDA00034498590100002213
Calculating the initial gaze tilt rate +.>
Figure BDA0003449859010000231
And a first corrected gaze inclination rate of change +>
Figure BDA0003449859010000232
According to the formula
Figure BDA0003449859010000233
Calculating a first corrected line of sight tilt angle theta Lp,i
If the current segment is the last segment of the straight line segment, L is i+1 0, thus the first corrected gaze inclination change rate
Figure BDA0003449859010000234
Directly taking 0.
Line-of-sight inclination rate of a current segment start point
Figure BDA0003449859010000235
And a first corrected gaze inclination rate of change +>
Figure BDA0003449859010000236
Averaging to obtain the average inclination angle change rate +.>
Figure BDA0003449859010000237
Whereby according to the initial time of flight prediction value +.>
Figure BDA0003449859010000238
Correcting the sight inclination angle of the current segment to obtain a first corrected sight inclination angle theta Lp,i
Correcting the sight line deflection angle of the current segment according to the following formula to obtain a first corrected sight line deflection angle phi Lp,i
If the current segment is in the turning segment, according to
Figure BDA0003449859010000239
Figure BDA00034498590100002310
Calculating the initial line of sight deflection change rate +.>
Figure BDA00034498590100002311
And a first corrected line-of-sight angle change rate
Figure BDA00034498590100002312
If the current segment is in the straight line segment, then according to +.>
Figure BDA00034498590100002313
Figure BDA00034498590100002314
Calculating the initial line of sight deflection change rate +.>
Figure BDA00034498590100002317
And a first corrected gaze angle change rate +>
Figure BDA00034498590100002315
According to the formula
Figure BDA00034498590100002316
Obtaining a first corrected sight offset angle psi Lp,i
If the current segment is the last segment of the straight line segment, L is i+1 Is 0, thus the first corrected line of sight offset angle change rate
Figure BDA0003449859010000241
Directly taking 0.
S2422, if the current segment is in the turn segment, according to equation y p,i =-D tc,i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, then according to formula y p,i =-L i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i
S2423, based on the first corrected gaze inclination angle θ Lp,i First corrected line-of-sight offset angle ψ Lp,i And a first speed predictive value v p,i Calculating the components of the first speed vector in all coordinate axes in the base reference coordinate system; calculating a first corrected ballistic inclination angle theta based on the components of each coordinate axis of the first velocity vector in the base reference coordinate system p,i
Specifically, the ballistic tilt angle θ is calculated according to the formula in step S231 as θ=f θLLyz V) correcting the inclination angle theta of the sight line Lp,i First corrected line-of-sight offset angle ψ Lp,i First speed predictive value v p,i Pretilt η of the end point of the current segment z,i+1 And the pre-bias angle eta of the end point of the current segment y,i+1 Carry in f θ Calculating to obtain a first corrected trajectory inclination angle theta p,i
Wherein the pretilt η of the end point of the current segment z,i+1 And the pre-bias angle eta of the end point of the current segment y,i+1 Can be calculated according to step S201.
S2424, based on the first speed predictive value v p,i First corrected ground height y p,i First corrected ballistic tilt angle theta p,i Calculating a first corrected speed change rate according to a kinetic equation of the mass center motion of the aircraft
Figure BDA0003449859010000242
Specifically, if the current segment is in the turn segment, the method is according to the formula
Figure BDA0003449859010000243
Calculating components of acceleration in a ballistic coordinate system; if the current segment is in the straight line segment, then according to the formula +. >
Figure BDA0003449859010000251
The components of acceleration in the ballistic coordinate system are calculated.
It should be noted that, if the current segment is the last segment of the straight line segment, then due to L i+1 0, so a is not calculated using the above formula yp,i And a zp,i Directly convert a yp,i And a zp,i Set to 0.
Wherein the pretilt η of the end point of the current segment z,i+1 Front offset angle eta of current segment end point y,i+1 Target distance D of current segment end point tc,i+1 Remaining straight-line segment voyage L of current segment end point i Can be calculated according to step S201.
According to the formula for calculating the speed change rate in steps S233-S237
Figure BDA0003449859010000252
First corrected ballistic inclination angle theta p,i First speed predictive value v p,i 、a yp,i 、a zp,i And a first corrected ground height y p,i Carry in->
Figure BDA0003449859010000253
Calculating to obtain a first corrected speed change rate->
Figure BDA0003449859010000254
S243, according to the formula
Figure BDA0003449859010000255
Calculating to obtain a first average correction speed +.>
Figure BDA0003449859010000256
/>
I.e. for initial rate of speed change
Figure BDA0003449859010000257
And a first corrected speed change rate->
Figure BDA0003449859010000258
Averaging to obtain an average speed change rate
Figure BDA0003449859010000259
According to the initial velocity v i Average rate of speed change->
Figure BDA00034498590100002510
And an initial time-of-flight prediction value->
Figure BDA00034498590100002511
Correcting the speed to obtain a first average corrected speed +.>
Figure BDA00034498590100002512
The accuracy of the estimation is improved by taking the average value.
S244, if the current segment is in the turning segment, according to
Figure BDA00034498590100002513
Calculating the current segment time of flight correction +. >
Figure BDA00034498590100002514
If the current segment is in the straight line segment, then according to +.>
Figure BDA00034498590100002515
Calculating the current segment time of flight correction +.>
Figure BDA00034498590100002516
According to the first average correction speed
Figure BDA0003449859010000261
Recalculating the time of flight of the current segment to obtain a current time of flight correction value +.>
Figure BDA0003449859010000262
The predicted value of the flight time is further corrected by correcting the line-of-sight inclination angle, the line-of-sight deflection angle and the trajectory inclination angle according to the initial speed change rate, so that the prediction is more accurate.
Step S25, based on the state variable of the current segment start point and the current segment flight time correction value
Figure BDA0003449859010000263
Predicting a state variable for a current segment end point, comprising:
s251, according to the initial speed change rate of the current segment
Figure BDA0003449859010000264
And said current segment time of flight correction value +.>
Figure BDA0003449859010000265
According to the formula->
Figure BDA0003449859010000266
Obtaining a second speed predicted value v of the current segment end point q,i
And calculating the speed predicted value again by using the corrected flight time estimated value, thereby improving the estimation accuracy.
S252, second speed predicted value v based on the segment end point q,i And the current segment time of flight correction value
Figure BDA0003449859010000267
Calculating a second rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure BDA0003449859010000268
In particular, according to
Figure BDA0003449859010000269
Calculating a second corrected line of sight tilt angle θ Lq,i According to
Figure BDA00034498590100002610
Calculating to obtain a second sight offset angle phi Lq,i The method comprises the steps of carrying out a first treatment on the surface of the According to y q,i =-D tc,i+1 sinθ Lq,i Calculating a second corrected ground height y q,i
The ballistic tilt angle θ is calculated according to the formula in step S231 as θ=f θLLyz V), the second velocity prediction value v q,i Second corrected line of sight tilt angle theta Lq,i Second corrected line-of-sight offset angle ψ Lq,i Pretilt η of the end point of the current segment z,i+1 And the pre-bias angle eta of the end point of the current segment y,i+1 Carry in f θ Calculating to obtain a second corrected trajectory inclination angle theta q,i
According to calculation a in step S2424 yp,i And a zp,i Will be v p,i Replaced by v q,i Calculating to obtain a yq,i And a zq,i
According to the formula for calculating the speed change rate in steps S233-S237
Figure BDA0003449859010000271
Second corrected ballistic tilt angle theta q,i Second speed predictor v q,i 、a yq,i 、a zq,i And a second corrected ground height y q,i Carry in->
Figure BDA0003449859010000272
Calculating a second corrected speed change rate +.>
Figure BDA0003449859010000273
S253, for initial speed change Rate
Figure BDA0003449859010000274
And a second corrected speed change rate->
Figure BDA0003449859010000275
Averaging to obtain a second average correction speed +.>
Figure BDA0003449859010000276
Figure BDA0003449859010000277
/>
S254, based on the second average correction speed
Figure BDA0003449859010000278
Calculating a third speed predicted value of the end point of the current segment
Figure BDA0003449859010000279
The third speed predicted value is the predicted speed v of the ending point of the current segment i+1
S255, according to the sight inclination angle theta of the current segmentation starting point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value
Figure BDA00034498590100002710
Calculating the sight inclination angle theta of the end point of the current segment L,i+1 The method comprises the steps of carrying out a first treatment on the surface of the Line-of-sight offset angle psi according to the current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value +.>
Figure BDA00034498590100002711
Calculating to obtain the sight offset angle psi of the end point of the current segment L,i+1
Specifically, the line-of-sight inclination angle θ of the end point of the current segment is calculated according to the following formula L,i+1
If the current segment is in the turning segment, according to the formula
Figure BDA00034498590100002712
Calculating a second corrected gaze angle change rate +.>
Figure BDA00034498590100002713
If the current segment is in the straight line segment, according to the formula +.>
Figure BDA00034498590100002714
Calculating a second corrected gaze angle change rate +.>
Figure BDA00034498590100002715
According to the formula
Figure BDA00034498590100002716
Calculating the line-of-sight inclination angle theta of the end point of the current segment L,i+1
If the current segment is in the turning segment, according to the formula
Figure BDA0003449859010000281
Calculating a second corrected gaze angle change rate +.>
Figure BDA0003449859010000282
If the current segment is in the straight line segment, according to the formula +.>
Figure BDA0003449859010000283
Calculating a second corrected gaze angle change rate +.>
Figure BDA0003449859010000284
According to the formula
Figure BDA0003449859010000285
Calculating the line-of-sight offset angle psi of the end point of the current segment L,i+1
S256, if the current segment is in the turning segment, according to the formula y i+1 =-D tc,i+1 sinθ L,i+1 Calculating the ground height of the ending point of the current segment; if the current segment is in a straight lineSegment according to y i+1 =-L i+1 sinθ L,i+1 Calculating the ground height y of the end point of the current segment i+1
By adopting the segmentation mode, the flight time is predicted for each segment, and the time sum of all segments is the total residual flight time, so that the estimation of the residual flight time is more accurate under the conditions that the speed of the unpowered aircraft is obviously reduced and the total leading angle is greatly changed.
The result of the remaining time of flight estimation using the prediction method and other methods of the present application is shown in fig. 4. It can be seen that the actual remaining time of flight is a straight line with a slope of-1, and the conventional estimation method and the constant-speed small-angle estimation method obviously estimate a large deviation because the speed variation and the large lead angle are not considered. By adopting the residual flight time prediction method, the residual flight time estimation curve is approximately overlapped with the actual residual flight time curve, the approach degree is higher and higher along with the time, and higher estimation precision is displayed.
In another aspect, an embodiment of the present invention provides a system for predicting a remaining time of flight of an unpowered aircraft, including:
the data acquisition module is used for acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of the aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight line flight of the current flight behavior according to the total lead angle;
the turning flight prediction module is used for segmenting the turning section according to the total front angle if turning flight is performed, and predicting the remaining flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
And the linear flight prediction module is used for calculating a remaining linear segment course according to the total lead angle and the target distance if the linear flight is linear flight, segmenting the linear segment according to the remaining linear segment course, and predicting the remaining flight time of the linear segment by adopting segmentation iteration, wherein the remaining flight time of the linear segment is the predicted total remaining flight time.
The method embodiment and the system embodiment are based on the same principle, and the related parts can be mutually referred to and can achieve the same technical effect. The specific implementation process refers to the foregoing embodiment, and will not be described herein.
In the invention, the technical schemes can be mutually combined to realize more preferable combination schemes. Additional features and advantages of the invention will be set forth in the description which follows, and in part will be obvious from the description, or may be learned by practice of the invention. The objectives and other advantages of the invention may be realized and attained by the structure particularly pointed out in the written description and drawings.
Those skilled in the art will appreciate that all or part of the flow of the methods of the embodiments described above may be accomplished by way of a computer program to instruct associated hardware, where the program may be stored on a computer readable storage medium. Wherein the computer readable storage medium is a magnetic disk, an optical disk, a read-only memory or a random access memory, etc.
The present invention is not limited to the above-mentioned embodiments, and any changes or substitutions that can be easily understood by those skilled in the art within the technical scope of the present invention are intended to be included in the scope of the present invention.

Claims (4)

1. A method of predicting the remaining time of flight of an unpowered aircraft, comprising the steps of:
acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of an aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight flight of the current flight behavior according to the total lead angle;
if the turning flight is the turning flight, segmenting the turning section according to the total lead angle, and predicting the residual flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
If the straight line flight is performed, calculating a remaining straight line segment course according to the total lead angle and the target distance, segmenting the straight line segment according to the remaining straight line segment course, and predicting the remaining flight time of the straight line segment by adopting segmented iteration, wherein the remaining flight time of the straight line segment is the predicted total remaining flight time;
predicting the residual flight time of a turning section by adopting piecewise iteration or predicting the residual flight time of a straight line section by adopting piecewise iteration comprises the following steps:
for each segment, a state variable of the starting point of the current segment is obtained, wherein the state variable comprises a velocity module value v i Inclination angle theta of sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i Ground level y i
If the current segment is in the turning segment, the state variable further comprises a total lead angle eta i And target distance D tc,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, the state variables further comprise the remaining straight line segment voyage L i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000011
Calculating the initial flight time predicted value of the current segment +.>
Figure FDA0004231250740000012
If the current segment is in the straight line segment, then according to +.>
Figure FDA0004231250740000013
Calculating the current segmentStart time of flight prediction +.>
Figure FDA0004231250740000014
Wherein Δη i Representing the variation of the total lead angle of the current segment, delta L representing the range variation of the remaining straight line segment of the current segment, K N Representing a navigation ratio;
calculating the initial speed change rate of the current segment according to the dynamics equation of the mass center motion of the aircraft based on the state variable of the starting point of the current segment
Figure FDA0004231250740000021
Based on the initial rate of speed change
Figure FDA0004231250740000022
An initial time-of-flight prediction value for said current segment>
Figure FDA0004231250740000023
Correcting to obtain a current segment flight time correction value +.>
Figure FDA0004231250740000024
State variable based on current segment start point and said current segment time-of-flight correction value
Figure FDA0004231250740000025
Predicting a state variable of a current segment end point, wherein the state variable of the current segment end point is a state variable of a next segment start point;
all segmented time of flight correction values
Figure FDA0004231250740000026
And remaining time of flight as a predicted turn or straight segment;
based on the initial rate of speed change
Figure FDA0004231250740000027
An initial time-of-flight prediction value for said current segment>
Figure FDA0004231250740000028
Correcting to obtain a current segment flight time correction value +.>
Figure FDA0004231250740000029
Comprising the following steps:
based on the initial rate of speed change
Figure FDA00042312507400000210
And said initial time of flight prediction value +.>
Figure FDA00042312507400000211
According to the formula->
Figure FDA00042312507400000212
Obtaining a first speed predictive value v of the segment end point p,i
A first speed predictor v based on the segment end point p,i And the initial time-of-flight prediction value
Figure FDA00042312507400000213
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid >
Figure FDA00042312507400000214
According to the formula
Figure FDA00042312507400000215
Calculating to obtain a first average correction speed +.>
Figure FDA00042312507400000216
If the current segment is in the turning segment, according to
Figure FDA00042312507400000217
Calculating the current segment time of flight correction +.>
Figure FDA00042312507400000218
If the current segment is in the straight line segment, then according to +.>
Figure FDA00042312507400000219
Calculating a current segment time-of-flight correction value
Figure FDA00042312507400000220
A first speed predictor v based on the segment end point p,i And the initial time-of-flight prediction value
Figure FDA0004231250740000031
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure FDA0004231250740000032
Comprising the following steps:
according to the first speed predictive value v p,i And the initial time-of-flight prediction value
Figure FDA0004231250740000033
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i
If the current segment is in the turning segment, then according to formula y p,i =-D tc,i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, then according to formula y p,i =-L i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the Wherein D is tc, i +1 Representing segment i+1Target distance of starting point, L i+1 Representing the rest straight-line segment course of the starting point of the i+1 segment;
based on the first corrected line of sight inclination angle theta Lp,i First corrected line-of-sight offset angle ψ Lp,i And a first speed predictive value v p,i Calculating the components of each coordinate axis of the first speed vector in the base reference coordinate system; calculating a first corrected ballistic inclination angle theta based on the components of each coordinate axis of the first velocity vector in the base reference coordinate system p,i
Based on the first speed predictive value v p,i First corrected ground height y p,i First corrected ballistic tilt angle theta p,i Calculating a first speed change rate according to a kinetic equation of the mass center motion of the aircraft
Figure FDA0004231250740000034
According to the first speed predictive value v p,i And the initial time-of-flight prediction value
Figure FDA0004231250740000035
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i Comprising:
correcting the sight inclination angle of the current segment according to the following formula to obtain a first corrected sight inclination angle theta Lp,i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000041
Calculating the initial gaze tilt rate +.>
Figure FDA0004231250740000042
And a first corrected gaze inclination rate of change +>
Figure FDA0004231250740000043
If the current segment is in the straight line segment, then according to
Figure FDA0004231250740000044
Calculating the initial gaze tilt rate +.>
Figure FDA0004231250740000045
And a first corrected gaze inclination rate of change +>
Figure FDA0004231250740000046
According to the formula
Figure FDA0004231250740000047
Obtaining a first corrected sight inclination angle theta Lp,i
Correcting the sight line deflection angle of the current segment according to the following formula to obtain a first corrected sight line deflection angle phi Lp,i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000048
Figure FDA0004231250740000049
Calculating the initial line of sight deflection change rate +.>
Figure FDA00042312507400000410
And a first corrected line-of-sight angle change rate
Figure FDA00042312507400000411
If the current segment is in the straight line segment, then according to
Figure FDA00042312507400000412
Figure FDA00042312507400000413
Calculating the initial line of sight deflection change rate +. >
Figure FDA00042312507400000414
And a first corrected gaze angle change rate +>
Figure FDA00042312507400000415
According to the formula
Figure FDA00042312507400000416
Obtaining a first corrected sight offset angle psi Lp,i
State variable based on current segment start point and said current segment time-of-flight correction value
Figure FDA0004231250740000051
Predicting a state variable for a current segment end point, comprising:
initial rate of speed change according to current segment
Figure FDA0004231250740000052
And said current segment time of flight correction value +.>
Figure FDA0004231250740000053
According to the formula
Figure FDA0004231250740000054
Obtaining a second speed predicted value v of the current segment end point q,i
A second velocity predictor v based on the segment end point q,i And the current segment time of flight correction value
Figure FDA0004231250740000055
Calculating a second rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure FDA0004231250740000056
For initial rate of speed change
Figure FDA0004231250740000057
And second rate of change of speed->
Figure FDA0004231250740000058
Averaging to obtain a second average correction speed +.>
Figure FDA0004231250740000059
Based on the second average corrected speed
Figure FDA00042312507400000510
Calculating a third speed predicted value of the end point of the current segment
Figure FDA00042312507400000511
The third speed predicted value is the predicted speed v of the ending point of the current segment i+1
Line of sight inclination angle theta according to current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value
Figure FDA00042312507400000512
Calculating the sight inclination angle theta of the end point of the current segment L,i+1 The method comprises the steps of carrying out a first treatment on the surface of the Line-of-sight offset angle psi according to the current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value +.>
Figure FDA00042312507400000513
Calculating to obtain the sight offset angle psi of the end point of the current segment L,i+1
If the current segment is in the turning segment, then according to formula y i+1 =-D tc,i+1 sinθ L,i+1 Calculating the current segmentThe ground height of the end point; if the current segment is in a straight line segment, then according to y i+1 =-L i+1 sinθ L,i+1 Calculating the ground height of the ending point of the current segment;
calculating the initial speed change rate of the current segment according to the dynamics equation of the mass center motion of the aircraft based on the state variable of the starting point of the current segment
Figure FDA00042312507400000514
Comprising the following steps:
based on the inclination angle theta of the sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i And velocity module v i Calculating the components of the speed vector in each coordinate axis in the base reference coordinate system; calculating the component of each coordinate axis in the basic reference coordinate system based on the velocity vector to obtain the initial ballistic inclination angle theta of the current segment i
If the current segment is in the turning segment, the method is based on the formula
Figure FDA0004231250740000061
Figure FDA0004231250740000062
Calculating components of acceleration in a ballistic coordinate system; if the current segment is in the straight line segment, the method is according to the formula
Figure FDA0004231250740000063
Calculating components of acceleration in a ballistic coordinate system;
according to the ground height y i Interpolation to obtain the atmospheric density ρ and the sound speed aT, and Mach number is calculated according to the atmospheric density ρ and the sound speed aT
Figure FDA0004231250740000064
And dynamic pressure->
Figure FDA0004231250740000065
Calculating the lift coefficient C according to the following formula Ld,i
F y =MT·g tc cosθ i +MT·a y,i
F z =MT·a z,i
Figure FDA0004231250740000066
Figure FDA0004231250740000067
Based on Mach and lift coefficient C Ld,i Calculating the attack angle command alpha by Newton iteration method c,i
Based on the Mach number Mach and the angle of attack command alpha c,i Calculating a resistance coefficient C according to the fitting formula D,i
According to the formula
Figure FDA0004231250740000068
The rate of change of the initial velocity is calculated,
wherein MT represents the mass of the aircraft, ST represents the reference area of the aircraft, g tc Representing the gravitational constant.
2. The method of claim 1, wherein the initial rate of change of velocity is based on a state variable of a current segment start point obtained for each segment
Figure FDA0004231250740000071
An initial time-of-flight prediction value for said current segment>
Figure FDA0004231250740000072
Before the correction, the method further comprises the following steps:
if the current segment is in the turning segment, thenCalculating the pretilt angle eta of the current segment end point according to the following formula z,i+1 Front deflection angle eta y,i+1 And a total lead angle eta i+1
Figure FDA0004231250740000073
Figure FDA0004231250740000074
Figure FDA0004231250740000075
η y,i+1 =η y,i +Δη y,i
η z,i+1 =η z,i +Δη z,i
η i+1 =arccos(cosη y,i+1 cosη z,i+1 );
According to the formula
Figure FDA0004231250740000076
Calculating a target distance D of the end point of the current segment tc,i+1
If the current segment is in the straight line segment, calculating the remaining straight line segment range L of the ending point of the current segment according to the following formula i+1 Front tilt angle eta z,i+1 And a prepositive deflection angle eta y,i+1
L i+1 =L i -ΔL
Figure FDA0004231250740000077
Figure FDA0004231250740000078
Wherein K is N Represents the navigation ratio, delta η Indicating the pre-angle segment quantity constant and NumP indicates the number of segments of the turn segment.
3. The method of claim 1, wherein calculating a total lead angle based on the line-of-sight inclination, line-of-sight offset, and velocity vector comprises:
calculating a direction cosine matrix from the reference coordinate system to the sight line coordinate system according to the sight line inclination angle and the sight line deflection angle;
calculating a component of the velocity vector in the line-of-sight coordinate system based on the direction cosine matrix;
the total lead angle is calculated according to the following formula:
Figure FDA0004231250740000081
Figure FDA0004231250740000082
η=arccos(cosη y cosη z )
wherein v is xL 、v yL And v zL Respectively representing the components of the velocity vector in three coordinate axes of a sight line coordinate system, eta z Represents the prepositive dip angle eta y Represents the lead angle, and η represents the total lead angle.
4. A system for predicting the remaining time of flight of an unpowered aircraft, comprising the following modules:
the data acquisition module is used for acquiring a target distance, a sight inclination angle, a sight deflection angle and a speed vector of the aircraft at the current moment, calculating a total lead angle based on the sight inclination angle, the sight deflection angle and the speed vector, and judging turning flight or straight line flight of the current flight behavior according to the total lead angle;
the turning flight prediction module is used for segmenting the turning section according to the total front angle if turning flight is performed, and predicting the remaining flight time of the turning section by adopting segmented iteration; calculating the remaining linear segment voyage according to the total lead angle of the starting point of the linear segment and the target distance, segmenting the linear segment according to the remaining linear segment voyage, and predicting the remaining flight time of the linear segment by adopting segmented iteration; the sum of the residual flight time of the turning section and the residual flight time of the straight line section is the predicted total residual flight time;
The linear flight prediction module is used for calculating a remaining linear segment course according to the total lead angle and the target distance if the linear flight is performed, segmenting the linear segment according to the remaining linear segment course, and predicting the remaining flight time of the linear segment by adopting segmented iteration, wherein the remaining flight time of the linear segment is the predicted total remaining flight time;
the turning flight prediction module adopts piecewise iteration to predict the remaining flight time of a turning section or adopts piecewise iteration to predict the remaining flight time of a straight line section, and comprises the following steps:
for each segment, a state variable of the starting point of the current segment is obtained, wherein the state variable comprises a velocity module value v i Inclination angle theta of sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i Ground level y i
If the current segment is in the turning segment, the state variable further comprises a total lead angle eta i And target distance D tc,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, the state variables further comprise the remaining straight line segment voyage L i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000091
Calculating the initial flight time predicted value of the current segment +.>
Figure FDA0004231250740000092
If the current segment is in the straight line segment, then according to +.>
Figure FDA0004231250740000093
Calculating the initial flight time predicted value of the current segment +.>
Figure FDA0004231250740000094
Wherein Δη i Representing the variation of the total lead angle of the current segment, delta L representing the range variation of the remaining straight line segment of the current segment, K N Representing a navigation ratio;
calculating the initial speed change rate of the current segment according to the dynamics equation of the mass center motion of the aircraft based on the state variable of the starting point of the current segment
Figure FDA0004231250740000095
Based on the initial rate of speed change
Figure FDA0004231250740000096
An initial time-of-flight prediction value for said current segment>
Figure FDA0004231250740000097
Correcting to obtain a current segment flight time correction value +.>
Figure FDA0004231250740000098
State variable based on current segment start point and said current segment time-of-flight correction value
Figure FDA0004231250740000099
Predicting a state variable of a current segment end point, wherein the state variable of the current segment end point is a state variable of a next segment start point;
all segmented time of flight correction values
Figure FDA0004231250740000101
And remaining time of flight as a predicted turn or straight segment;
based on the initial rate of speed change
Figure FDA0004231250740000102
An initial time-of-flight prediction value for said current segment>
Figure FDA0004231250740000103
Correcting to obtain a current segment flight time correction value +.>
Figure FDA0004231250740000104
Comprising the following steps:
based on the initial rate of speed change
Figure FDA0004231250740000105
And said initial time of flight prediction value +.>
Figure FDA0004231250740000106
According to the formula->
Figure FDA0004231250740000107
Obtaining a first speed predictive value v of the segment end point p,i
A first speed predictor v based on the segment end point p,i And the initial time-of-flight prediction value
Figure FDA0004231250740000108
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid >
Figure FDA0004231250740000109
According to the formula
Figure FDA00042312507400001010
Calculating to obtain a first average correction speed +.>
Figure FDA00042312507400001011
If the current segment is in the turning segment, according to
Figure FDA00042312507400001012
Calculating the current segment time of flight correction +.>
Figure FDA00042312507400001013
If the current segment is in the straight line segment, then according to +.>
Figure FDA00042312507400001014
Calculating a current segment time-of-flight correction value
Figure FDA00042312507400001015
A first speed predictor v based on the segment end point p,i And the initial time-of-flight prediction value
Figure FDA00042312507400001016
Calculating a first rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure FDA00042312507400001017
Comprising the following steps:
according to the first speed predictive value v p,i And the initial time-of-flight prediction value
Figure FDA00042312507400001018
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i
If the current segment is in the turning segment, then according to formula y p,i =-D tc,i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the If the current segment is in a straight line segment, then according to formula y p,i =-L i+1 sin(θ Lp,i ) Calculating a first corrected ground height y p,i The method comprises the steps of carrying out a first treatment on the surface of the Wherein D is tc,i+1 Representing the target distance of the starting point of the i+1 segment, L i+1 Representing the rest straight-line segment course of the starting point of the i+1 segment;
based on the first corrected line of sight inclination angle theta Lp,i First corrected line-of-sight offset angle ψ Lp,i And a first speed predictive value v p,i Calculating the components of each coordinate axis of the first speed vector in the base reference coordinate system; calculating a first corrected ballistic inclination angle theta based on the components of each coordinate axis of the first velocity vector in the base reference coordinate system p,i
Based on the first speed predictive value v p,i First corrected ground height y p,i First corrected ballistic tilt angle theta p,i Calculating a first speed change rate according to a kinetic equation of the mass center motion of the aircraft
Figure FDA0004231250740000111
According to the first speed predictive value v p,i And the initial time-of-flight prediction value
Figure FDA0004231250740000112
Correcting the sight inclination angle and the sight deflection angle of the current segment to respectively obtain a first corrected sight inclination angle theta Lp,i And a first corrected line-of-sight offset angle ψ Lp,i Comprising:
correcting the sight inclination angle of the current segment according to the following formula to obtain a first corrected sight inclination angle theta Lp,i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000113
Calculating the initial gaze tilt rate +.>
Figure FDA0004231250740000114
And a first corrected gaze inclination rate of change +>
Figure FDA0004231250740000115
If the current segment is in the straight line segment, then according to
Figure FDA0004231250740000116
Calculating the initial gaze tilt rate +.>
Figure FDA0004231250740000117
And a first corrected gaze inclination rate of change +>
Figure FDA0004231250740000118
According to the formula
Figure FDA0004231250740000119
Obtaining a first corrected sight inclination angle theta Lp,i
Correcting the sight line deflection angle of the current segment according to the following formula to obtain a first corrected sight line deflection angle phi Lp,i
If the current segment is in the turning segment, according to
Figure FDA0004231250740000121
Figure FDA0004231250740000122
Calculating the initial line of sight deflection change rate +.>
Figure FDA0004231250740000123
And a first corrected gaze angle change rate +>
Figure FDA0004231250740000124
If the current segment is in the straight line segment, then according to
Figure FDA0004231250740000125
Figure FDA0004231250740000126
Calculating the initial line of sight deflection change rate +. >
Figure FDA0004231250740000127
And a first corrected gaze angle change rate +>
Figure FDA0004231250740000128
According to the formula
Figure FDA0004231250740000129
Obtaining a first corrected sight offset angle psi Lp,i
State variable based on current segment start point and said current segment time-of-flight correction value
Figure FDA00042312507400001210
Predicting a state variable for a current segment end point, comprising:
initial rate of speed change according to current segment
Figure FDA00042312507400001211
And said current segment time of flight correction value +.>
Figure FDA00042312507400001212
According to the formula
Figure FDA00042312507400001213
Obtaining a second speed predicted value v of the current segment end point q,i
A second velocity predictor v based on the segment end point q,i And the current segment time of flight correction value
Figure FDA00042312507400001214
Calculating a second rate of change of velocity from the kinetic equation of the movement of the aircraft centroid>
Figure FDA00042312507400001215
For initial rate of speed change
Figure FDA00042312507400001216
And second rate of change of speed->
Figure FDA00042312507400001217
Averaging to obtain a second average correction speed +.>
Figure FDA00042312507400001218
Based on the second average corrected speed
Figure FDA0004231250740000131
Calculating a third speed predicted value of the end point of the current segment
Figure FDA0004231250740000132
The third speed predicted value is the predicted speed v of the ending point of the current segment i+1
Line of sight inclination angle theta according to current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value
Figure FDA0004231250740000133
Calculating the sight inclination angle theta of the end point of the current segment L,i+1 The method comprises the steps of carrying out a first treatment on the surface of the Line-of-sight offset angle psi according to the current segment start point L,i Speed v of end point of current segment i+1 And said current segment time of flight correction value +.>
Figure FDA0004231250740000134
Calculating to obtain the sight offset angle psi of the end point of the current segment L,i+1
If the current segment is in the turning segment, then according to formula y i+1 =-D tc,i+1 sinθ L,i+1 Calculating the ground height of the ending point of the current segment; if the current segment is in a straight line segment, then according to y i+1 =-L i+1 sinθ L,i+1 Calculating the ground height of the ending point of the current segment;
calculating the initial speed change rate of the current segment according to the dynamics equation of the mass center motion of the aircraft based on the state variable of the starting point of the current segment
Figure FDA0004231250740000135
Comprising the following steps:
based on the inclination angle theta of the sight line L,i Sight line deflection angle psi L,i Front tilt angle eta z,i Front deflection angle eta y,i And velocity module v i Calculating the components of the speed vector in each coordinate axis in the base reference coordinate system; calculating the component of each coordinate axis in the basic reference coordinate system based on the velocity vector to obtain the initial ballistic inclination angle theta of the current segment i
If the current segment is in the turning segment, the method is based on the formula
Figure FDA0004231250740000136
Figure FDA0004231250740000137
Calculating components of acceleration in a ballistic coordinate system; if the current segment is in the straight line segment, the method is according to the formula
Figure FDA0004231250740000138
Calculating components of acceleration in a ballistic coordinate system;
according to the ground height y i Interpolation to obtain the atmospheric density ρ and the sound speed aT, and Mach number is calculated according to the atmospheric density ρ and the sound speed aT
Figure FDA0004231250740000141
And dynamic pressure->
Figure FDA0004231250740000142
Calculating the lift coefficient C according to the following formula Ld,i
F y =MT·g tc cosθ i +MT·a y,i
F z =MT·a z,i
Figure FDA0004231250740000143
Figure FDA0004231250740000144
Based on Mach and lift coefficient C Ld,i Calculating the attack angle command alpha by Newton iteration method c,i
Based on the Mach number Mach and the angle of attack command alpha c,i Calculating a resistance coefficient C according to the fitting formula D,i
According to the formula
Figure FDA0004231250740000145
The rate of change of the initial velocity is calculated,
wherein MT represents the mass of the aircraft, ST represents the reference area of the aircraft, g tc Representing the gravitational constant.
CN202111672205.0A 2021-12-31 2021-12-31 Method and system for predicting residual flight time of unpowered aircraft Active CN114326813B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111672205.0A CN114326813B (en) 2021-12-31 2021-12-31 Method and system for predicting residual flight time of unpowered aircraft

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111672205.0A CN114326813B (en) 2021-12-31 2021-12-31 Method and system for predicting residual flight time of unpowered aircraft

Publications (2)

Publication Number Publication Date
CN114326813A CN114326813A (en) 2022-04-12
CN114326813B true CN114326813B (en) 2023-06-20

Family

ID=81020045

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111672205.0A Active CN114326813B (en) 2021-12-31 2021-12-31 Method and system for predicting residual flight time of unpowered aircraft

Country Status (1)

Country Link
CN (1) CN114326813B (en)

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104077469A (en) * 2014-05-28 2014-10-01 中国人民解放军海军航空工程学院 Speed prediction based segmentation iteration remaining time estimation method
CN107132765A (en) * 2017-06-01 2017-09-05 烟台南山学院 A kind of angle-of-attack based on trajectory planning and attack time control method
CN109484674A (en) * 2018-10-12 2019-03-19 湖北航天技术研究院总体设计所 A kind of real-time track maneuver autopilot method based on target track parameter
CN113671974A (en) * 2021-07-18 2021-11-19 北京理工大学 Turning approach accurate guidance method for return section of cross-domain aircraft

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR2942566B1 (en) * 2009-02-24 2016-01-22 Thales Sa METHOD FOR MANAGING THE FLIGHT OF AN AIRCRAFT

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104077469A (en) * 2014-05-28 2014-10-01 中国人民解放军海军航空工程学院 Speed prediction based segmentation iteration remaining time estimation method
CN107132765A (en) * 2017-06-01 2017-09-05 烟台南山学院 A kind of angle-of-attack based on trajectory planning and attack time control method
CN109484674A (en) * 2018-10-12 2019-03-19 湖北航天技术研究院总体设计所 A kind of real-time track maneuver autopilot method based on target track parameter
CN113671974A (en) * 2021-07-18 2021-11-19 北京理工大学 Turning approach accurate guidance method for return section of cross-domain aircraft

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Impact time control using biased proportional navigation for missiles with varying velocity;Guoxin SUN 等;Chinese Journal of Aeronautics;全文 *
基于速度预测的导引律剩余时间估计;张友安 等;北京航空航天大学学报;第43卷(第10期);全文 *

Also Published As

Publication number Publication date
CN114326813A (en) 2022-04-12

Similar Documents

Publication Publication Date Title
CN114326814B (en) Three-dimensional guidance system of unpowered aircraft
US10656650B2 (en) Method for guiding and controlling drone using information for controlling camera of drone
US20120123615A1 (en) Method and a system for estimating a trajectory of a moving body
WO2020000127A1 (en) Navigation path tracking control method, device, mobile robot and system
CN112498744B (en) Longitudinal and transverse loose coupling online track planning method and electronic equipment
CN106406333B (en) A kind of stratospheric airship pitch angle tracking based on integral form terminal sliding mode
CN110702106A (en) Unmanned aerial vehicle, course alignment method and device thereof and storage medium
CN106840194A (en) A kind of Large azimuth angle linear alignment method
CN104121930B (en) A kind of compensation method based on the MEMS gyro drift error adding table coupling
CN114326813B (en) Method and system for predicting residual flight time of unpowered aircraft
CN102607557B (en) GPS/IMU (Global Position System/Inertial Measurement Unit)-based direct integral correction method for aircraft attitudes
CN114020019A (en) Guidance method and device for aircraft
CN117389312A (en) Model-based three-dimensional tracking control method for counter roll of underwater vehicle
CN116301058B (en) Unmanned flight feedback nonlinear yaw control method, system and equipment
CN113110527A (en) Cascade control method for finite time path tracking of autonomous underwater vehicle
CN116796119A (en) Monocular distance measurement precision calculation method based on unmanned aerial vehicle motion platform
CN112034879A (en) Standard trajectory tracking guidance method based on height-range ratio
CN110647161A (en) Under-actuated UUV horizontal plane trajectory tracking control method based on state prediction compensation
CN115268475A (en) Robot fish accurate terrain tracking control method based on finite time disturbance observer
CN114660587A (en) Jump and glide trajectory target tracking method and system based on Jerk model
CN102508819B (en) Angular-speed-based quaternion Legendre approximate output method during extreme flying of aircraft
CN102495830B (en) Quaternion Hartley approximate output method based on angular velocities for aircraft during extreme flight
CN102495831B (en) Quaternion Hermitian approximate output method based on angular velocities for aircraft during extreme flight
KR102114051B1 (en) Method of non-linear control for aircraft considering center of gravity movement
CN102495829B (en) Quaternion Walsh approximate output method based on angular velocities for aircraft during extreme flight

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