Skip Re-Entry Trajectory Detection and Guidance for Maneuvering Vehicles

The re-entry trajectory of maneuvering vehicles with medium to high hypersonic lift-to-drag ratios is generally planned using quasi-equilibrium flight conditions known from Space Shuttles. They may exhibit an oscillation re-entry phenomenon termed skip re-entry when related components or sensors fail. However, conventional re-entry guidance only considers quasi-equilibrium flights and ignores the possibility of the occurrence of an unexpected skip trajectory; this may lead to the failure of the re-entry mission due to a lack of a corresponding guidance strategy. However, the detection of a skip trajectory is the necessary reference for the decision-making of calling a related guidance algorithm that helps improve the safety of vehicle re-entry. Herein, a skip re-entry detection and trajectory control solution is proposed to play an emergency role in the cases of skip re-entry. Firstly, the oscillation frequency characteristics of the linearized re-entry motion equation of a vehicle are analyzed, and an approximate analytical relationship is constructed for skip altitude estimation. Then, the residual deviation between the altitude feedback data and the estimated skip altitude is calculated and compared with the threshold to determine the occurrence of skip re-entry. In addition, a method for controlling the skip re-entry trajectory with the range extension is developed by controlling the bank angle with a fixed angle of attack profile, satisfying the path constraint requirements. The results indicate that the method effectively performs skip re-entry detection and that it can help extend the range of the vehicles in abnormal re-entry scenarios, keeping the flight within the path constraints and guiding it to the expected location.


Introduction
Experiences of the Space Shuttle provide a baseline for planning and optimizing the re-entry trajectory for vehicles with similar symmetric configurations and medium or high lift-to-drag ratios [1,2]. From low Earth orbit (LEO) back into the atmosphere, the drag-acceleration profile planning algorithm [3] is explicit and practical. Several research studies [4][5][6][7][8][9] are also based on the planning of drag-acceleration profiles, for which the acceleration sensor plays an indispensable role. However, electromagnetic radiation in the near-earth space environment is very complex [10] and may cause sensor performance degradation or failure [11]. However, sensor performance degradation is not easy to detect directly and will lead to a deviation of the result of the guidance algorithm calculation, which may lead to an abnormal re-entry trajectory.
For the improved safety and probability of success, different deviations and abnormal scenarios need to be taken into consideration during vehicle development. For commercial space launch and re-entry activities, the Federal Aviation Administration requires that the trajectory analysis for both normal and malfunction flights must characterize variability and include deviation in paragraphs 450. 117  of maneuvering vehicles with medium to high hypersonic lift-to-drag ratios is generally planned by employing quasi-equilibrium flight conditions known from the Space Shuttle [3]. In addition, there are different planning techniques to adapt to various scenarios. Shen et al. provided a rapid generation method of three degrees of freedom (3-DOF) re-entry [12]. Subsequently, Lu provided a unified re-entry guidance method, independent of the hypersonic lift-to-drag ratio (L/D) of the vehicle, which can range from 0.28 for a capsule crew exploration vehicle to 3.5 for a common aero vehicle [13], basically all flying in quasi-equilibrium glide mode. Another possibility, however, is to fly along non-equilibrium trajectories, such as skips that may occur under abnormal conditions. The skip re-entry trajectory planning was first introduced for the Apollo lunar mission to achieve a long downrange from the orbits; this was done to provide diversion capability to avoid bad and disruptive weather conditions [14]. Orion is also designed with the capability of performing a skip entry for extending the down range [15]. For the skip entry guidance theory, Loh [16,17] proposed a first and second order method to provide an analytical approximate solution, and then Vinh et al. [18] refined and further developed the second order solution for lifting skip trajectories.
For re-entry vehicles, other skipping strategies, such as reference-following controllers [19,20] and numeric predictor-corrector skip re-entry algorithms have been proposed to provide onboard real-time trajectory generation and guidance [21,22]. Furthermore, Brunner and Lu [23,24] integrated trajectory prediction with closed-loop correcting guidance using the bank angle determined based on the downrange requirement. Subsequently, Luo et al. [25,26] presented a skip guidance algorithm that used a numerical predictor-corrector and a patched corridor for low-lifting capsules returning from the Moon, where the distribution of the bank angle was piecewise and linear with respect to the normalized energy. Cheng et al. [27] proposed a numerical multi-constrained predictor-corrector guidance algorithm, which designed a bank corridor to help convert the trajectory planning into a root-finding problem and developed a constraint management module to improve the satisfaction of path and terminal constraints. Liu et al. [28] represented an attempt to apply a second order cone programming, a branch of convex optimization, to a class of highly nonlinear trajectory optimization problem in entry flight. In addition, Wang et al. [29] developed the convex optimization methods to solve hypersonic trajectory optimization problems, which generated a reference trajectory by solving a second order cone programming problem and designed an optimal feedback guidance law using a constrained quadratic programming method to track the trajectory. Moreover, skipping guidance is also closely related to aerocapture, such as the analytic predictor-corrector guidance, which is being discussed for the Mars Sample Return Orbiter vehicle [30], the numerical optimal predictor-corrector guidance algorithms that are being developed to guide the spacecraft through aerocapture into a target orbit [31,32].
For re-entry of a maneuvering vehicle, if only conventional re-entry trajectory guidance is considered and the possibility of skip trajectory is ignored, there may be a loss of control during re-entry due to a lack of any corresponding guidance strategy. During the process of LEO re-entry, sensor data plays a key role in the guidance algorithm, and its failure or performance degradation may lead to a skip trajectory. On one hand, during the process of initial re-entry, the sensor data provides a judgment basis for using different guidance algorithms [13,14,23,24]. When a failure or the performance degradation of any sensor leads to a system misjudgment, the guidance instructions will follow the initial descent instructions. Then, the vehicle goes to a skip mode due to its natural phugoid oscillatory tendency in a dense atmosphere [14] and experiences serious deviation from the original trajectory planning. On the other hand, a quasi-equilibrium glide trajectory planning strategy is widely used for re-entry [12,13]. In this method, the bank angle is normally controlled to suppress the occurrence of a skip; however, any sensor performance degradation may lead to a guidance output which does not restrain the skip, and this may lead to a failure of the guidance algorithm in the skip flight.
It can be seen that, when there is a skip in a routine re-entry process, calling on the corresponding skip trajectory guidance strategy will help with improving the safety of the re-entry vehicle, and the Sensors 2020, 20, 2976 3 of 25 detection of a skip trajectory will be a necessary reference for decision-making when calling the guidance algorithm.
However, conventional sensors used in skip re-entry detection provide only basic data and fail to further determine whether the vehicle will undergo skip re-entry. In addition, for the above skip guidance algorithms, which are mainly employed for expected skip guidance scenarios, such as aerocapture [30][31][32], and extended range entry [15], there is no skip trajectory detection and transition logic, and they are mainly focused on vehicles with L/D ratios below 0.5 at re-entry and velocities close to the escape velocity of the Earth, that is, 11.2 km/s [19][20][21][22][23][24][25][26]30]. However, for the LEO re-entry, a conventional and practical approach for re-entry is to glide in a quasi-equilibrium condition. On this basis, this paper considers an unexpected skip occurrence and expects to monitor or assist decision-making through skip trajectory detection and calls on the corresponding skip guidance algorithm to improve the safety of LEO re-entry.
Therefore, this article proposes a skip re-entry detection and trajectory control solution to deal with the skip problem in trajectory control during the normal re-entry of hypersonic maneuvering winged vehicles with a medium or high L/D.
First, the oscillation frequency characteristics of the linearized equation of motion of the vehicle re-entry are analyzed, and an approximate relationship between the spatial frequency and altitude is obtained. Second, the residual deviation between the altitude feedback data and the estimated skip altitude is calculated and compared with a threshold value to determine whether the vehicle will undergo skip re-entry. In addition, a skip re-entry trajectory control method with the range extension, which is realized by controlling the bank angle of a profile with a fixed angle of attack, and a numerical search algorithm are discussed. Finally, the performance is assessed by skip re-entry detection rate tests, false alarm rate tests and Monte Carlo dispersion simulations.

Equations and Constraints of Re-Entry Flight
For re-entry, the 3-DOF equations, including the Earth's curvature and rotation [33], are widely accepted, as follows: where r is the radial distance from the center of Earth to the vehicle, normalized by the radius of the Earth r e = 6378 km, θ and φ are the longitude and latitude, respectively, γ is the flight path angle relative to the surface of the Earth, ψ is the course angle describing the relative velocity vector measured clockwise from the north, σ is the bank angle measured positive to the right from the view inside the vehicle, V is the velocity relative to Earth, m is the mass of the vehicle, g is the gravitational acceleration, Ω e is the rotation rate of Earth, L and D are aerodynamic lift and drag, respectively, as expressed below: where ρ is the local atmospheric density used in this study, which is in accordance with the 1976 U.S. standard atmosphere model [34], C L and C D are lift and drag coefficients, respectively, S re f is the reference area of the vehicle. In re-entry, the path constraints often include limits of heat rate, normal load factor, dynamic pressure, typically described as [13] .
Q is the heat rate at a stagnation point on the surface of the vehicle, .
Q max is the heat rate limit, and C Q is a constant parameter. N Z and n max are the normal load factor and its limit, and q and q max are the dynamic pressure and its maximum value.
If in a quasi-equilibrium glide condition (QEGC), the parameters retain the following relation [12,35]: where σ QEGC is the bank angle. This constraint reduces oscillations and preserves the bank angle margin. This is a soft constraint in the sense that its enforcement need not be too strict.
The re-entry vehicle is expected to reach a desired terminal condition, for which the state is specified by the altitude h * f and velocity V * f at a distance R * f from the target location. The final state of re-entry is described by the final altitude h f and velocity V f at a distance R f from the target location, and the terminal constraints are given by The exponential atmosphere model [34] used for constraint transformation is: where h = r − r e , h r is the atmospheric density scale height and ρ 0 = 1.225 kg/m 3 . The tracking is realized when the velocity vs. altitude plane in Equations (3)-(6) is transformed by replacing the density with the altitude in Equation (8) into By analyzing the influence of the angle of attack (AoA) on the heat rate constraint, normal load factor constraint and dynamic pressure constraint, the shape of the AoA profile can be designed. In addition, as a result of the success of the path constraint management in the Shuttle strategy, an AoA profile can be pre-planned based on the nominal range-to-go, maximum L/D, etc. [3], with an example, as follows: (13) where the maximum AoA α max is selected to avoid over-limit heat rate in rapid deceleration at initial re-entry, α L/Dmax is the AoA at maximum L/D for final glide, V α_max and V α_min are the maximum and minimum switch velocities, V is the current velocity, k represents the parameter of gradient between α L/D max and α max , and a denotes a compensation constant.

Skip Re-Entry Detection Solution
For low hypersonic L/D vehicles whose re-entry velocity is close to the escape velocity of 11.2 km/s, the Kepler phase often exists because of the high re-entry energy state. Conversely, high hypersonic L/D vehicles may not enter the Kepler phase without sufficient re-entry energy from low orbits. In this study, it is assumed that the skip re-entry of the vehicle is carried out in the atmosphere, and the detection of skip is focused on longitudinal motion.

Analysis on Phugoid Oscillation
As the equilibrium glide motion is performed at a baseline where . γ ≈ 0 and the rotation rate of the Earth is neglected, the longitudinal state vector x = [r V γ] T and the control vector u = cos σ are linearized [36] according to Equations (1), (2) and (9), as follows: where D r , D V , L r , and L V are the partial derivatives. γ assumes small values in the re-entry process of equilibrium glide for medium and high L/D vehicles [37]. Then, the characteristic equations are the following: According to Etkin [38] and Laitone and Chou [39], there are three modes for longitudinal motion at hypersonic speed, namely, the phugoid mode, short-period mode and spiral mode, and the root of the spiral mode is real and close to zero. Therefore, the "phugoid oscillation" mode is estimated by ignoring the C 3 term, as follows: Then the "phugoid oscillation" frequency is approximated as and its corresponding spatial frequency is During the skips occurring at hypersonic Mach numbers, the AoA is given from the profile at α = α max , and the lift and drag coefficients are almost constant; therefore, the phugoid frequency changes with flight altitude and velocity. Figure 1 presents a typical frequency variation in a skipping re-entry. From the figure, it can be seen that the frequency is comparable to the ripple of the trajectory; furthermore, the frequency increases as the velocity decreases and shows additional oscillation with oscillations in the altitude. In view that the change in the frequency reflects the basic shape of the altitude vs. the velocity, an approximate skip altitude h skip is constructed with the piecewise frequency, as follows: where h max is the predicted peak skipping altitude, and F is a constant factor determined according to the starting skip state, as follows: where h 0 and Ω 0 are the starting altitude and the corresponding spatial frequency, respectively.
and its corresponding spatial frequency is During the skips occurring at hypersonic Mach numbers, the AoA is given from the profile at max = , and the lift and drag coefficients are almost constant; therefore, the phugoid frequency changes with flight altitude and velocity. Figure 1 presents a typical frequency variation in a skipping re-entry. From the figure, it can be seen that the frequency is comparable to the ripple of the trajectory; furthermore, the frequency increases as the velocity decreases and shows additional oscillation with oscillations in the altitude. In view that the change in the frequency reflects the basic shape of the altitude vs. the velocity, an approximate skip altitude skip h is constructed with the piecewise frequency, as follows: where max h is the predicted peak skipping altitude, and F is a constant factor determined according to the starting skip state, as follows: where 0 h and 0  are the starting altitude and the corresponding spatial frequency, respectively.   (21) and (22) and numerical integration (the Runge-Kutta method is employed) as a baseline from Equations (1), (2), and (13). The construction function approximation matches the numerical integration, while the deviation of the construction function approximation accumulates with continuous skips. To reduce the deviation, the initial state parameters 0 h , 0  and max h in Equation (22) can then be updated at each skip. Figure 2b,c show that the accumulated deviation is significantly reduced for the second and third skip, if the initial state parameters are updated for these skips.   (21) and (22) and numerical integration (the Runge-Kutta method is employed) as a baseline from Equations (1), (2), and (13). The construction function approximation matches the numerical integration, while the deviation of the construction function approximation accumulates with continuous skips. To reduce the deviation, the initial state parameters h 0 , Ω 0 and h max in Equation (22) can then be updated at each skip. Figure 2b,c show that the accumulated deviation is significantly reduced for the second and third skip, if the initial state parameters are updated for these skips.

Skip Detection Based on Estimated Skip Altitude
The residual deviation between the estimated skip altitude (h skip ) and the altitude output of the sensor measurement feedback (h) is used as the observation result in each period to detect the skip, that is: Skip condition: A means test is used for the residual evaluation for a selected skip velocity scope, and the analysis is performed in Section 5.2.1.
Establishing binary hypothesis: H 0 : Skip Define the false alarm rate P F = P(H 1 / H 0 ), missing detection rate P M = P(H 0 / H 1 ), and detection rate P D = 1 − P M . The designed skip re-entry detection function is: where w k is the weight coefficient and A k is the h k group mean. According to the central limit theorem, the mean value follows the normal distribution, and then the weighted data of the mean value follows the standard normal distribution, λ k ∼ N(0, 1). The skip is defined as The no-skip is defined as To determine T D according to the requirement of the system P F : then, where α T is the expected value.

Skip Re-Entry Detection for Trajectory Control Logic
For vehicles with high hypersonic L/D, the tendency for "phugoid oscillation" in re-entry trajectory has been suppressed by augmenting guidance [13]. However, the oscillations have a range-increasing capability and can be used to extend the downrange.
The range rate is From Equations (1) and (2), Equation (29) is written as The hypersonic aerodynamic coefficients are almost constant with the Mach number at a specified AoA. As the flight path tends to be flat, e.g., in an altitude change of 40 km compared with a range of 2000 km, it is assumed that γ ≈ 0 for the range estimation. To compare the range extending capability Sensors 2020, 20, 2976 9 of 25 between quasi-equilibrium glide flights and skip flights, it is assumed that they are at the same latitude φ = 0 • . Then, Equation (30) can be rewritten as where V 1 and V 2 are the specified velocities for concerned re-entry, normally It is shown that for a given vehicle and scheduled AoA profile, the atmospheric density is a main factor affecting re-entry range. In addition, by suppressing skips, cos σ QEGC < cos σ skip , compared with the quasi-equilibrium glide flight, the skip motion at the same velocity will fly at a thinner atmospheric density, and thus ρ QEGC > ρ skip . The quasi-equilibrium glide detours continuously by alternative banking and corresponds to a lower downrange than the skip: After detecting the skip re-entry, the characteristics of the skip with increasing range can be used to guide the vehicle to the desired alternate landing area, and the skip range can be predicted by numerical integration from Equations (1), (2), (13), and (29). The transition logic from skip re-entry detection to trajectory control is shown in Figure 3. The hypersonic aerodynamic coefficients are almost constant with the Mach number at a specified AoA. As the flight path tends to be flat, e.g., in an altitude change of 40 km compared with a range of 2000 km, it is assumed that 0   for the range estimation. To compare the range extending capability between quasi-equilibrium glide flights and skip flights, it is assumed that they are at the same latitude = 0


. Then, Equation (30) can be rewritten as where 1 V and 2 V are the specified velocities for concerned re-entry, normally α_max It is shown that for a given vehicle and scheduled AoA profile, the atmospheric density is a main factor affecting re-entry range. In addition, by suppressing skips, QEGC cos cos skip   , compared with the quasi-equilibrium glide flight, the skip motion at the same velocity will fly at a thinner atmospheric density, and thus QEGC skip   . The quasi-equilibrium glide detours continuously by alternative banking and corresponds to a lower downrange than the skip: After detecting the skip re-entry, the characteristics of the skip with increasing range can be used to guide the vehicle to the desired alternate landing area, and the skip range can be predicted by numerical integration from Equations (1), (2), (13), and (29). The transition logic from skip re-entry detection to trajectory control is shown in Figure 3.

Skip Re-Entry Trajectory Control
In this section, the implementation of the trajectory control by generating suitable bank angles when skip occurs is performed. The constraints are the heat rate, normal load factor and the dynamic pressure considered.
The re-entry is normally divided into the down control phase and final glide phase [12]; however, once the skip occurs, the skip trajectory needs to be controlled. Therefore, in this study, the skip re-entry is divided into the down control, skip control, and the final glide phases, in a similar manner to the Apollo re-entry [14], as shown in Figure 4.

Skip Re-Entry Trajectory Control
In this section, the implementation of the trajectory control by generating suitable bank angles when skip occurs is performed. The constraints are the heat rate, normal load factor and the dynamic pressure considered.
The re-entry is normally divided into the down control phase and final glide phase [12]; however, once the skip occurs, the skip trajectory needs to be controlled. Therefore, in this study, the skip re-entry is divided into the down control, skip control, and the final glide phases, in a similar manner to the Apollo re-entry [14], as shown in Figure 4. The down control phase is the same as in normal re-entry where the vehicle flies through vacuum into the atmosphere from the initial re-entry interface to the state when the descent rate is close to zero. During the down control phase, the bank angle down  is normally given by a constant which needs to satisfy boundary constraints [3].
The different bank angle control methods are used in skip control phase and final glide phase.
With the  versus flight velocity profile given, a constant bank angle skip  is employed in the skip control phase and reference trajectory-tracking method is used in the final glide. For the skip re-entry, the range threshold specified energy state is a reference condition, and the control phase transition logic is used for bank angle control method selection in the skip control phase and the final glide phase.
In normal re-entry, it is assumed that the quasi-equilibrium glide flight is employed for guidance, which suppress the skips. Therefore, we can obtain a skip guidance command magnitude constraint with the quasi-equilibrium glide flight condition [12], as follows: In addition, when the AoA profile and the magnitude of constant down  are preset, the skip trajectory needs to satisfy the constraint boundary Equations (9)- (11). To predict the skip motion from a velocity-altitude profile, the Runge-Kutta numerical iteration method using Equations (1), (2), (13) and (9)  Method 1: The reference trajectory-tracking algorithm To further satisfy the path constraints, a 3rd order reference baseline trajectory is designed for the final glide phase, as follows: The down control phase is the same as in normal re-entry where the vehicle flies through vacuum into the atmosphere from the initial re-entry interface to the state when the descent rate is close to zero. During the down control phase, the bank angle σ down is normally given by a constant which needs to satisfy boundary constraints [3].
The different bank angle control methods are used in skip control phase and final glide phase. With the α versus flight velocity profile given, a constant bank angle σ skip is employed in the skip control phase and reference trajectory-tracking method is used in the final glide. For the skip re-entry, the range threshold specified energy state is a reference condition, and the control phase transition logic is used for bank angle control method selection in the skip control phase and the final glide phase.
In normal re-entry, it is assumed that the quasi-equilibrium glide flight is employed for guidance, which suppress the skips. Therefore, we can obtain a skip guidance command magnitude constraint with the quasi-equilibrium glide flight condition [12], as follows: In addition, when the AoA profile and the magnitude of constant σ down are preset, the skip trajectory needs to satisfy the constraint boundary Equations (9)- (11). To predict the skip motion from a velocity-altitude profile, the Runge-Kutta numerical iteration method using Equations (1), (2), (13) and (9)-(11) can be employed to estimate the next skip state for a given σ skip . Then, the maximum skip guidance command magnitude σ skip_max can be further identified within the dispersion tests of atmospheric density, lift and drag coefficients, and mass. Based on Equation (33), the skip guidance command magnitude constraint can be obtained as follows: Method 1: The reference trajectory-tracking algorithm To further satisfy the path constraints, a 3rd order reference baseline trajectory is designed for the final glide phase, as follows: where k i (i = 0, 1, 2, 3) are solved using reference points P1, P2, P3, and P4. P1 is the final point of the down control phase and P4 is the terminal constraint point, including the final velocity and altitude conditions. P2 is the state point with zero descent rate after the first skip that can be predicted based on the Runge-Kutta method. P3 is a designed point between P2 and P4 to construct the reference baseline trajectory. All the points are needed to satisfy the velocity-altitude path constraints profile, as shown in Figure 5. The skip altitude for the skip motion is higher than that of the reference baseline trajectory at the same velocity; this further proves that the skip satisfies path constraints.
Sensors 2020, 20, x FOR PEER REVIEW 11 of 27 where ( 0,1, 2,3) i ki  are solved using reference points P1, P2, P3, and P4. P1 is the final point of the down control phase and P4 is the terminal constraint point, including the final velocity and altitude conditions. P2 is the state point with zero descent rate after the first skip that can be predicted based on the Runge-Kutta method. P3 is a designed point between P2 and P4 to construct the reference baseline trajectory. All the points are needed to satisfy the velocity-altitude path constraints profile, as shown in Figure 5. The skip altitude for the skip motion is higher than that of the reference baseline trajectory at the same velocity; this further proves that the skip satisfies path constraints. To track the reference baseline trajectory and obtain the behavior of a stable second order feedback system in altitude deviation, Equation (1) is rewritten as: and, From the above description, the velocity vs. altitude profile can be used as a monitor. Then, the minimum magnitude constraint of guidance command can be obtained when Equation (6)   To track the reference baseline trajectory and obtain the behavior of a stable second order feedback system in altitude deviation, Equation (1) is rewritten as: and, ..
Because h re f belongs to the velocity-altitude profile Equations (9) h, and the altitude deviation feedback control law can be designed as the following, where the parameter ω h and ξ h are the designed undamped natural frequency and damping. Inserting Equations (1), (36), and (37) into Equation (38) leads to the following expression: From the above description, the velocity vs. altitude profile can be used as a monitor. Then, the minimum magnitude constraint of guidance command can be obtained when Equation (6) takes the minimum value by σ EQ , and the maximum magnitude constraint needs to satisfy the path constraints. Hence, we have where σ path is the bank angle according to the path constraints and can be obtained using the path constraints velocity and altitude in equilibrium glide condition, where . γ ≈ 0 and the rotation rate of the Earth is neglected: Therefore, if σ f inal is beyond the boundary, a preset boundary value satisfying the Equation (40) will be employed for guidance command.

Method 2: The control phase transition logic
Based on the range-to-go, the skip trajectory is adjusted through control phase transition logic, including the selection of tracking h re f and employing skips by the threshold range R thres , as follows: and, The range threshold is the boundary to decide the transition from the skip control phase to the final control phase depending on the skip ending velocity V 2 , as shown in Figure 6. The interval from V 1 to V 2 is employed for "phugoid oscillation", and the range boundary R thres is approximately converted into the accumulated range after tracking h re f from V 2 to the end. V 1 is the velocity at the state when the descent rate is close to zero, V α_max is a designed velocity in Equation (13). The range threshold state is specified as the related altitude h thres and velocity V thres at a distance R thres from the target location.
The range threshold is the boundary to decide the transition from the skip control phase to the final control phase depending on the skip ending velocity   For the lateral guidance logic, the sign of the bank angle is determined by the threshold for bank angle reverse, which is used during the flight phases. In the skip control phase, the reverse threshold is the deviation of the cross range. In the final glide phase, the deviations of the cross range and course angle are employed to determine the sign of the bank angle for different distance conditions, as follows:

Method 3: The lateral guidance logic
For the lateral guidance logic, the sign of the bank angle is determined by the threshold for bank angle reverse, which is used during the flight phases. In the skip control phase, the reverse threshold is the deviation of the cross range. In the final glide phase, the deviations of the cross range and course angle are employed to determine the sign of the bank angle for different distance conditions, as follows: |∆c r |≤ ∆c thres R togo ≥ R course ∆ψ ≤ ∆ψ thres R togo < R course (44)  where ∆c r and ∆ψ are the deviations of the cross range and course angle, respectively. ∆c thres and ∆ψ thres are the respective deviation thresholds of the cross range and course angle, ∆c thres is designed based on the velocity function, and the course angle is controlled by entering the range boundary R course . R togo is the distance between the current position of the vehicle and the destination position. When the deviation is higher than the threshold, the bank angle sign is changed. The ∆ψ thres selecting principle is mainly based on the terminal course angle constraint, and the ∆c thres choosing rule is determined based on the bank angle reversal frequency. In this paper, σ skip is a main variable and is used for adjusting the extension range in the re-entry process. To establish a skip re-entry trajectory control, the range threshold specified states, V * thres , h * thres , and e * thres = 1 2 V * 2 thres + gh * thres , are preset in a way to satisfy the final re-entry condition (the footprint method can be used to design and evaluate the specified state [40]), and then σ skip is calculated to guide the vehicle to the desired range threshold state.
The σ skip numerical search algorithm is presented in Algorithm 1. for (i = 0 to N) do 6: 3-DOF equations numerical integration output assignment 7: Compare R thres and R togo , call Method 2; 8: if R thres < R togo then 9: α = α pro f il , σ =σ skip , call Method 3; 10: else 11: call Method 1; 12: To form skip re-entry guidance, α and σ are the control variables in each guidance cycle, while α is provided by the profile represented by Equation (13). The calculation of σ depends on different guidance strategies of Equations (39), (42), (43) and (44) and the σ skip numerical search algorithm for each flight phase. During the down control phase, σ skip can be predicted to prepare for skip motion.
When a skip re-entry is detected, the above numerical search algorithm can be used to correct σ skip according to the current state in a designed period. At the final glide phase, the reference trajectory-tracking method is employed to enhance the robustness.

Example Vehicle and Preset Parameters
The example vehicle is taken to be similar to the Space Shuttle [41], with a mass of 76,297 kg and a hypersonic lift-to-drag ratio of 1.884 at a Mach number of 20 and an angle of attack of 40 • . The pre-planned AoA profile corresponds to Equation (13), as follows: (V ≥ 5000 m/s) 0.00625V + 15 (1000 m/s ≤ V < 5000 m/s) 15 (V < 1000 m/s) When the flight vehicle returns to the atmosphere from a low Earth orbit, the deorbit impulse and transfer orbit are mainly calculated according to the given re-entry window into the atmosphere (including initial re-entry altitude, velocity, and flight path angle). Therefore, the assumed initial longitude and latitude are given in this paper, the initial re-entry altitude, velocity, and flight path angle are provided according to references [42,43], and the mission conditions are listed in Table 1.

Mean Test
In this study, residual deviations between the estimated skip altitude and the altitude output of the sensor measurement feedback are compared in each period and used as a skip re-entry detection condition. Figure 7b,d show the mean residual deviations between the estimated skip altitude and the sensor feedback data. Because the altitude of the skip is generally higher than that of quasi-equilibrium glide, Figure 7b shows a small mean residual deviation when a skip re-entry happens; however, Figure 7d shows a significant negative deviation when the re-entry is in quasi-equilibrium glide. Therefore, a mean threshold was used as the decision-making condition to unambiguously distinguish the skip re-entry state from the no skip re-entry state.
the sensor feedback data. Because the altitude of the skip is generally higher than that of quasi-equilibrium glide, Figure 7b shows a small mean residual deviation when a skip re-entry happens; however, Figure 7d shows a significant negative deviation when the re-entry is in quasi-equilibrium glide. Therefore, a mean threshold was used as the decision-making condition to unambiguously distinguish the skip re-entry state from the no skip re-entry state. (c) skip altitude estimation vs. the sensor altitude feedback data during quasi-equilibrium glide re-entry; (d) quasi-equilibrium glide re-entry residuals between the skip altitude estimation and the sensor altitude feedback data.

False Alarm Rate and Detection Rate Test
Based on the description in Section 3.2, the false alarm rate and the detection rate of the skip re-entry are defined as F P and D P . In this study, it is assumed that the process noise is zero mean white noise with a standard deviation of 1/2 Q , and the standard deviation of the noise changes from  (d) quasi-equilibrium glide re-entry residuals between the skip altitude estimation and the sensor altitude feedback data.

False Alarm Rate and Detection Rate Test
Based on the description in Section 3.2, the false alarm rate and the detection rate of the skip re-entry are defined as P F and P D . In this study, it is assumed that the process noise is zero mean white noise with a standard deviation of Q 1/2 , and the standard deviation of the noise changes from 1 × 10 −6 to 1 × 10 −2 . Each simulation generated 500 data points, and every five consecutive data points form a group to calculate the mean value. To design w k = 1/1800, a false alarm rate α T = 8%, and from the look-up table T D = −1.41 is employed. The test results for the false alarm rate are shown in Figure 8. Sensors 2020, 20, x FOR PEER REVIEW 16 of 27 Figure 8. False alarm rates test. Figure 8 shows that the false alarm rate varies from 3.5% to 6%, suggesting that the threshold determined by the formula can meet the requirement = 8% T  . Figure 9 shows that when the above values of the threshold and noise interference parameters are used in this detection rate test, the skip re-entry detection rate varies from 88% to 91.5%. Figure 9. Skip re-entry detection rates test. The computer configuration is an Intel(R) Core(TM) i5-4200H CPU @ 2.80GHz, RAM 4.00GB, Windows 8 operating system.  Figure 8 shows that the false alarm rate varies from 3.5% to 6%, suggesting that the threshold determined by the formula can meet the requirement α T = 8%. Figure 9 shows that when the above values of the threshold and noise interference parameters are used in this detection rate test, the skip re-entry detection rate varies from 88% to 91.5%.  Figure 8 shows that the false alarm rate varies from 3.5% to 6%, suggesting that the threshold determined by the formula can meet the requirement = 8%

Skip Re-Entry Trajectory Control Test
T  . Figure 9 shows that when the above values of the threshold and noise interference parameters are used in this detection rate test, the skip re-entry detection rate varies from 88% to 91.5%. Figure 9. Skip re-entry detection rates test. The computer configuration is an Intel(R) Core(TM) i5-4200H CPU @ 2.80GHz, RAM 4.00GB, Windows 8 operating system. Figure 9. Skip re-entry detection rates test.

Execution Time Performance Test of the σ skip Numerical Search Algorithm
This section is used to evaluate the running time for the σ skip numerical search algorithm from three aspects: each motion state calculation, eachσ skip given calculation from current state to the expected range threshold, and σ skip searching calculation, as shown in Figure 10a-e.
The computer configuration is an Intel(R) Core(TM) i5-4200H CPU @ 2.80GHz, RAM 4.00GB, Windows 8 operating system.  Figure 10a shows that the algorithm requires less than 35 ms to recur from the current state to the next state through a numerical integration of the 3-DOF equations using the RK4 method, which can meet the calculation requirements of the on-board computer within a period of 40 ms. Figure 10b,c are the calculation times of the algorithm cycle after the assumed bank angle command ˆs kip  is given at different positions. Figure 10b is the calculation from the re-entry position to the desired range threshold thres R , and Figure 10c is the calculation from the skip detection position to the expected range threshold thres R . By comparing the two figures, it can be seen that the farther the distance from the current position of the vehicle to thres R is, the longer the calculation time of the algorithm will be. Figure 10b shows that the calculation time of the algorithm is less than 4.5 s, while that of Figure 10c is less than 3 s. Figure 10d,e are the times for the algorithm to calculate the skip bank angle command based on the distance from the current position of the vehicle to the range threshold thres R . Figure 10d shows the calculation times for the current position of the vehicle at the down control phase of re-entry, and Figure 10e shows the calculation times for the current position at the skip control phase. The comparison between the two figures shows that the average calculation time of the former is about 10 s, while that of the latter is about 5 s, which shows that the algorithm's online calculation time is shorter when the vehicle approaches to the range threshold. This suggests that the algorithm supports the requirement that the online correction time decreases with decreasing range-to-go.

Skip Re-Entry Trajectory Control Test under Monte Carlo Method
The mission conditions are listed in Table 1. For the 1000 trials of the Monte Carlo analysis, Gaussian-shaped or uniform dispersions are assumed, as summarized in Table 2. The analytical atmospheric density dispersion with respect to the 1976 U.S. Standard Atmosphere is used to imitate atmospheric uncertainties [44].  Figure 10a shows that the algorithm requires less than 35 ms to recur from the current state to the next state through a numerical integration of the 3-DOF equations using the RK4 method, which can meet the calculation requirements of the on-board computer within a period of 40 ms. Figure 10b,c are the calculation times of the algorithm cycle after the assumed bank angle commandσ skip is given at different positions. Figure 10b is the calculation from the re-entry position to the desired range threshold R thres , and Figure 10c is the calculation from the skip detection position to the expected range threshold R thres . By comparing the two figures, it can be seen that the farther the distance from the current position of the vehicle to R thres is, the longer the calculation time of the algorithm will be. Figure 10b shows that the calculation time of the algorithm is less than 4.5 s, while that of Figure 10c is less than 3 s. Figure 10d,e are the times for the algorithm to calculate the skip bank angle command based on the distance from the current position of the vehicle to the range threshold R thres . Figure 10d shows the calculation times for the current position of the vehicle at the down control phase of re-entry, and Figure 10e shows the calculation times for the current position at the skip control phase. The comparison between the two figures shows that the average calculation time of the former is about 10 s, while that of the latter is about 5 s, which shows that the algorithm's online calculation time is shorter when the vehicle approaches to the range threshold. This suggests that the algorithm supports the requirement that the online correction time decreases with decreasing range-to-go.

Skip Re-Entry Trajectory Control Test under Monte Carlo Method
The mission conditions are listed in Table 1. For the 1000 trials of the Monte Carlo analysis, Gaussian-shaped or uniform dispersions are assumed, as summarized in Table 2. The analytical atmospheric density dispersion with respect to the 1976 U.S. Standard Atmosphere is used to imitate atmospheric uncertainties [44].
This paper focuses on detecting the unexpected skip re-entry and guiding the vehicle to the planned recovery area in case of an emergency. The terminal conditions are as follows: If the ending position satisfies the velocity and altitude requirements and is within the precision circle of a radius of 30 km, the mission is considered to be a success; beyond the circle of radius of 50 km, it is considered to be a failure, and if in between, it is counted as a 50% success.

Parameters Distribution Max
Re-entry initial longitude, deg Gaussian 0.1 Re-entry initial latitude, deg Gaussian 0.1 Re-entry initial relative velocity, m/s Gaussian 5 Re-entry initial flight path angle, deg Gaussian 0.05 Re-entry initial heading angle, deg Gaussian 0.5 Atmospheric density Analytical 30% Mass, kg Uniform ±1% Figure 11 shows the histories under Monte Carlo dispersion (a subset of 1000 simulations). The skip re-entry detection threshold and parameter set is the same as in Section 5.2.2., the initial state of re-entry is given by the parameters of mission 1 in Table 1, the roll rate is less than 5/s, the range threshold is preset to R thres = 4000 km, the range threshold specified state is e * thres ( V * thres = 6846 m/s, h * thres = 78,000 m), and σ skip ∈ [−50, 50]. During the skip control phase, the bank angle correction period is 100 s.  The normal re-entry plan for a vehicle is to fly from the initial re-entry positon to the target site; when the skip re-entry is detected, the vehicle can be operated in accordance with the standby emergency plan to turn to the emergency site for recovery. In the down control phase, skip  numerical search algorithm is used to predict the bank angle command in the case of skip re-entry. In the skip control phase, the algorithm is used to correct the energy needed to reach the range threshold, and the algorithm call cycle is 100 s in these tests. There are two elements that are specifically dedicated to the robustness of the algorithm, the first of which is the correction of bank angle in the skip control phase to satisfy the range threshold expected energy, and the second one is involved in the reference trajectory tracking law in the final glide phase to restrain the influence of process disturbances. Figure 12a shows the terminal miss distances for all tests, in which the terminal points within a circle of a radius of 30 km are assumed to be successful in view of the next stage of terminal area guidance, and the points within the circles of radius 30-50 km are considered to be semi-successful. For the 1000 trial simulations with parameter dispersions, the success rate is above 90%. In addition, when the vehicle arrives at the final re-entry velocity condition, Figure 12b shows that the altitude deviation varies from 29.6 km to 31 km, which satisfies the 5%  requirement. The main purpose of this emergency solution is to recover the vehicle, and the velocity of the vehicle arriving at the position is 1000 m/s. Under this condition, the vehicle can hover and descend to help the ground find its position.  The normal re-entry plan for a vehicle is to fly from the initial re-entry positon to the target site; when the skip re-entry is detected, the vehicle can be operated in accordance with the standby emergency plan to turn to the emergency site for recovery. In the down control phase, σ skip numerical search algorithm is used to predict the bank angle command in the case of skip re-entry. In the skip control phase, the algorithm is used to correct the energy needed to reach the range threshold, and the algorithm call cycle is 100 s in these tests. There are two elements that are specifically dedicated to the robustness of the algorithm, the first of which is the correction of bank angle in the skip control phase to satisfy the range threshold expected energy, and the second one is involved in the reference trajectory tracking law in the final glide phase to restrain the influence of process disturbances. Figure 12a shows the terminal miss distances for all tests, in which the terminal points within a circle of a radius of 30 km are assumed to be successful in view of the next stage of terminal area guidance, and the points within the circles of radius 30-50 km are considered to be semi-successful. For the 1000 trial simulations with parameter dispersions, the success rate is above 90%. In addition, when the vehicle arrives at the final re-entry velocity condition, Figure 12b shows that the altitude deviation varies from 29.6 km to 31 km, which satisfies the ±5% requirement.

Skip Re-Entry Detection and Trajectory Control Application
The Proposed Solution in Abnormal Skip Re-Entry Emergency Scenarios Figure 13 shows the cases of emergency re-entry of a maneuvering vehicle along an abnormal skip trajectory. The vehicle carries out a planned QEGC normal re-entry from the sites P1 and P2, The main purpose of this emergency solution is to recover the vehicle, and the velocity of the vehicle arriving at the position is 1000 m/s. Under this condition, the vehicle can hover and descend to help the ground find its position. Figure 13 shows the cases of emergency re-entry of a maneuvering vehicle along an abnormal skip trajectory. The vehicle carries out a planned QEGC normal re-entry from the sites P1 and P2, respectively. However, due to a degradation of sensor performance, the vehicle does not enter the QEGC control mode due to the conditions such as "the total aerodynamic acceleration is greater than or equal to 1.52 m/s 2 " in reference [13] or "the load is greater than 0.05 g" in reference [23], but still flies according to the guidance commands of the capture atmosphere, and the vehicle flies along the abnormal skip trajectory. If the QEGC control algorithm is not called, there is a delay in calling it or is not suitable for skip trajectory control, the vehicle may be out of control, which threatens the safety of the ground. Possible points of failure are indicated by a red "x", as shown in Figure 13a. In order to reduce the risk after such a situation occurs, the proposed solution for skip re-entry detection and trajectory control can be used as a reference method to guide the vehicle to the emergency area, as shown in Figure 13a

The Application of QEGC Method and the Proposed Solution in Abnormal Skip Re-Entry Scenarios
This section is devoted to the analysis of the influence of the QEGC method and the method presented in this paper on the trajectory control of skip re-entry. Assuming that the performance of the sensor is degraded, QEGC control is called with a delay, the bank angle continues to use the Figure 13. Emergency re-entry of a maneuvering vehicle for normal and abnormal trajectory control: (a) lateral histories for normal and abnormal re-entry trajectories, (b) velocity-altitude histories for normal and abnormal re-entry trajectories.

The Application of QEGC Method and the Proposed Solution in Abnormal Skip Re-Entry Scenarios
This section is devoted to the analysis of the influence of the QEGC method and the method presented in this paper on the trajectory control of skip re-entry. Assuming that the performance of the sensor is degraded, QEGC control is called with a delay, the bank angle continues to use the down control phase command, as shown in Figure 14c, and thus there will be a skipping phenomenon. The method described in this paper, with the function of skip detection, can help to drive the vehicle to the desired alternate area. However, QEGC method is no longer suitable for skip re-entry control, which makes it impossible to smoothly control the vehicle to an alternate area, as shown in Figure 14a. Generally, QEGC re-entry is based on the designed nominal trajectory, and it contains no detection function for an unexpected skip trajectory. As shown in Figure 14b, in the process of hypersonic re-entry, when the skip occurs, the delayed call of the QEGC algorithm is likely to increase the short-term skip amplitude. Though the original QEGC algorithm emphasizes the Generally, QEGC re-entry is based on the designed nominal trajectory, and it contains no detection function for an unexpected skip trajectory. As shown in Figure 14b, in the process of hypersonic re-entry, when the skip occurs, the delayed call of the QEGC algorithm is likely to increase the short-term skip amplitude. Though the original QEGC algorithm emphasizes the suppression of a skip, and it may be successful in suppressing the skip within a certain period of time, the vehicle would by then have deviated too far from its nominal trajectory, which may lead to a failure of the re-entry mission.
The solution proposed in this paper adds to the detection function of skip re-entry, which can utilize the extended range abilities generated by a skip. When an unexpected skip occurs, it can actively use the skip ability to guide the vehicle to the desired recovery area and satisfy the path constraints.

Discussion
For the skip trajectory control, a range threshold R thres also can be found by referring to flight history data. In addition, to improve the accuracy of the final position for future work, the increase in the correction frequency of a given skip command and the improvement of the guidance algorithm in the final glide phase will be studied. The final re-entry condition in this paper refers to the energy management window condition. In future work, the question of how to deal with the area guidance problem of emergency alternate terminal will be further studied.

Conclusions
This paper is dedicated to the detection of unexpected skip re-entry and active use of the natural phugoid oscillations of medium or high L/D re-entry of the vehicle to guide it to the emergency area for recovery. The following conclusions can be drawn.
(1) An approximate analytical relationship is constructed for skip altitude estimation based on the oscillation frequency characteristic of the linearized re-entry motion equation of the vehicle. Based on the above analytical relationship, the skip re-entry detection method can be used as a standby tool of airborne monitoring in the form of software to alert the skipping during the re-entry process or to prepare to call other emergency trajectory control strategies; (2) Based on the pre-planned angle of attack profile, a phase separation of the down control, skip control, and final glide phases is employed in this paper. On this basis, a control phase transition logic-based on the range threshold under the velocity-altitude profile is proposed, which can smoothly connect the skip control phase and final glide phase, and provide support in calling related guidance algorithms at each phase to achieve a single or multiple skip re-entry to improve range capabilities; (3) Simulations further demonstrate that the proposed solution can achieve an expected detection rate, its running time is reasonable, and the trajectory control can satisfy path constraints and be robust under Monte Carlo dispersions. Finally, it has also been demonstrated that the method can guide the vehicle to an emergency area for recovery when it skips. Funding: This research received no external funding.