Identification Method of SUAV in Diving Phase Based on Flight Tests

. A mathematical model of the dive phase is an important research content for improving the accuracy of terminal control in the small unmanned aerial vehicle. The acquisition of the diving model poses new challenges, such as the small installation space, ultra-low ﬂying height of small suicide drones, short ﬂight time, strong coupling, less observable measurement, and elastic deformation of the wings during the drone dive phase. Based on the autoregressive moving average method, a multi-input multioutput noise term topology mathematical model is proposed in this paper. Through an improved least squares identiﬁcation method, the diving model in the ﬂight test is analyzed and veriﬁed. The identiﬁcation results of the diving model obtained by the proposed method are compared with the least squares method dive model. The results indicate that the mathematical model and identiﬁcation method proposed in this paper can eﬀectively obtain the parameters of the drone dive model.


Introduction
With the development of SUAV technology, the function of unmanned aerial vehicles (UAVs) has gradually expanded from a single function of scouting, interference, and damage assessment to the integration of scouting and strikes. During the execution of scouting and searching above a target area, the SUAV can immediately hit the target when a timesensitive target is found. Due to strong maneuverability and high attack precision (with a controllable range at the meter level), the SUAV can effectively destroy the target while carrying fewer combat damage units, thereby exhibiting high benefits in actual combat. Representative SUAVs are the Switchblade, the Coyote, and the Hero loitering munition.
During attacks, dive modeling is quite important in the end trajectory control. To enhance the controlling precision of SUAV during attacks, the diving model should reflect the SUAV's actual performance. Researchers mainly adopted wind tunnel tests and computational fluid dynamics (CFD) software for the establishment of mathematical models of conventional air vehicles [1]. Both methods require high professional skills and are time-consuming as well as costly. Currently, certain engineering calculation software of aerodynamic parameters exhibit poor computational accuracy, which can thereby affect the precision of the entire mathematical model. By utilizing the system response to specific input single utilization, the system identification method can rapidly and effectively acquire the mathematical model of the entire air vehicle, which has been extensively applied during flight tests of conventional air vehicles [2]. Since the introduction of "identification" by American scholar Zadeh in 1956 [3], system identification, a method of system modeling, has been continuously developed for over 60 years. Following booming, the system identification was rapidly applied to system modeling of air vehicles [4][5][6].
e aforementioned has become a mature modeling method in the conventional air vehicle field.
Recently, with the development of low-cost microelectromechanical system (MEMS) sensor technology, system identification has been utilized in extensive applications in the SUAV field. With the unscented Kalman filter (UKF), Mararu Naruoka et al. [7] performed system identification on a Mitsubishi MARS07AF.
is was an SUAV carrying MEMS devices in 2009. However, the identification results significantly differed from the wind tunnel test data. In 2010, Khaled S. Hatamleh et al. used an inertial measurement unit (IMU) on the UAV for the related model parameters estimation [8]. In 2011, Andrei Dorobantu et al. employed the frequency-domain identification technique for both longitudinal and transverse model identification of a low-sized fixedwing UAV [9]. In 2012, Aaron Wypyszynski conducted modelidentification flight tests on a small UAV and identified its aerodynamic parameters [10]. In 2012, Arnaud Koehl et al. established a 6-DOF model of GLMAV, a coaxial-rotor UAV, and performed parameter identification through Kalman filtering [11]. In 2015, Konstantin identified the four-rotor UAV model by using the recursive algorithm [12]. Kailong Liu et al. employed a recursive algorithm to investigate the battery temperature prediction method with regard to the battery safety problem [13,14]. In addition, with the improvement of hardware processing capabilities, methods such as machine learning and Kalman filtering are also used in model identification and parameter prediction of complex systems [15][16][17][18].
Model identification of small UAVs research is currently mainly focused on longitudinal or horizontal models in the cruise phase [19][20][21][22][23]. e research on the diving phase of SUAVs is mainly focused on control rate, guidance rate, and recognition of the diving model. In 2017, Zhou [24] conducted an integrated design for the guidance and controller of the suicide UAV. Moreover, the author verified the effectiveness of the system. He [25] proposed a hypersonic gliding vehicle spiral dive maneuver penetration concept and designed an adaptive proportional guidance rate for virtual portrait targets. Karpenko [26] investigated the evaluation method of attack aircraft's dive characteristics and analyzed the time and height loss of attack aircraft in different dive and pull modes. Wang [27] studied the guidance and control laws of antitank missiles with curved firing and self-seeking. Moreover, the author established a guidance and control law suitable for ballistic control. With the rapid expansion of military applications of suicide attacks by small UAVs and to improve the accuracy of UAV strikes, the problem of accurate modeling and identification of the diving model needs to be solved urgently. As an SUAV gradually approaches the target during diving, the position and attitude should be constantly corrected. is phase has strong nonlinear characteristics and is a multi-input multioutput system [8,9]. In general, due to installation space and landing weight limitations, testing equipment such as attack angle and sideslip angle sensors cannot be installed. Moreover, only a certain amount of data can be acquired. In addition, since the SUAV always flies at a low velocity, the flight is easily disturbed by many factors, such as gusts.
In this paper, a topological mathematical model of multiinput and multioutput noise terms for SUAV is first constructed. en, an improved multivariable augmented recursive least squares identification method is proposed, and the feasibility of the method is verified by the simulation model. Finally, the identification methods are compared and verified through flight tests of SUAV.

Mathematical Model
e diving phase of SUAV is a very complex and highly coupled nonlinear dynamic system. Its six-degree-of-freedom motion equations can be expressed as [28,29] where μ, ], and ω are body velocities along with the x-axis, yaxis, and z-axis. Parameters p, q, and r are angular rates along the x-axis, y-axis, and z-axis. Parameters ϕ, θ, and ψ are Euler angles along the x-axis, y-axis, and z-axis. Coefficients C X , C Y , C Z , C L , C M , and C N are aerodynamic force and moment coefficients. C b w is the coordinate transformation matrix, and m, S, and b are the quality, reference area, and wingspan of the small UAV, respectively. Parameter q is the air density, is the inertia matrix, with I xx being the rolling moment of inertia, I yy being the pitching moment of inertia, I zz being the yawing moment of inertia, and I xz and I zx are inertia products. e parameter V A is the airspeed, while μ A , ] A , and ω A are airspeed vector components along the x-axis, y-axis, and z-axis in the body frame. Parameter n is the rotational speed of the propeller, D is the diameter of the propeller, and C F T is the thrust coefficient. e aerodynamic force of SUAV has a very complicated nonlinear relationship with many factors such as the geometry, flight speed, aerodynamic angle, and even atmospheric temperature of the aircraft. During the flight test, factors such as aerodynamic angle and wing deformation do not have the measurement conditions. If all these factors are simultaneously taken into account, model identification and solution will be relatively difficult. Limited by the research funding, the abovementioned complex parameters can be hardly obtained. erefore, the ARMA model is used in this paper to solve this problem. is method treats the model as a gray box, which can fit the model according to the input and output data of the system. Moreover, it can predict the output features of the model with sufficient accuracy, thereby reducing the amount of calculation and the difficulty of the identification test.
e ARMA model is a common and high-precision time-series short-term prediction method. It exhibits extensive applications in many domains such as fault detection, controller model identification, sensor correction, and structural testing and identification [30][31][32][33]. e dynamic discrete system input u(k) and the output sampling value series y(k) can be expressed as the following linear difference equation: where a i (i � 1, 2, . . . , n a ) and b i (i � 0, 1, . . . , n b ) are the coefficients of autoregression items. e multi-input multioutput system of the SUAV, consisting of r-dimensional inputs and m-dimensional outputs, is displayed in Figure 1. e control input matrix and system output matrix are represented as According to the analysis of the characteristics during diving, the structure of the diving model can be expressed as where A 1 , A 2 , . . . , A n a denotes a 2 * 2 matrix to be identified, B 0 , B 1 , . . . , B n b is also a 2 * 2 matrix to be identified, and V(k) denotes a 2 * 2 colored noise matrix in the model indicating interference factors such as wind gusts and temperature during the test. Parameters U(k) and Y(k) are the control input matrix and the state output matrix, respectively. U(k) and Y(k) can be written as where ele(k) and ail(k) are the k-th PWM control input data in pitching and rolling channels, while q(k) and p(k) denote the k-th pitching angular rate and the rolling angular rate corresponding to control inputs, respectively.
To determine the structural parameters n a and n b , the one-step-ahead prediction value of the model output can be written as e parameter θ is calculated through least squares, while the output prediction error vector can be defined as e loss function of the model can be defined as 3. Algorithm e least squares method was proposed by Gauss to determine the orbit of Ceres. Within the method, the optimal functional matching was searched for by employing the sum minimization of squares of the errors. Currently, the least squares method is a common method for parameter estimation on the test data that exhibits extensive applications and popularization in many domains. Despite simplicity and practicability, the least squares method is no longer an unbiased estimation for the system including colored noise. erefore, given that an SUAV is a multi-input multioutput coupling system during the diving phase, while the colored noise imposed considerable disturbance, a multivariable recursive augmented least squares identification method is derived in this paper.

Standard Augmented Recursive Least Squares Method.
e standard recursive least squares identification algorithm can be expressed as During the actual small aircraft test process, due to the influence of gust and other environmental effects, the actual noise is not white noise. erefore, it is necessary to identify where ] ′ represents random noise with the nonzero mean, ] represents random noise with the mean value of 0, and c 1 , c 2 , ..., c n c represent the parameters of the noise model that require identification. e model can be expressed as Noise sequence ] ′ (k) can be replaced by the estimation of ](k) { }: where e noise sequence can be supplemented as the identified object. Parameters to be reformulated are as follows: 1 , a 2 , . . . , a n a , b 1 , b 2 , . . . , b n b , c 1 , c 2 , . . . , c n c T . (13) e model can be expressed as e calculation expression of the augmented least squares recursive method can be formulated as

Improved Multivariable Augmented Recursive Least Squares
Method. e multivariable identification model can be defined as where To improve the model's ability to deal with the noise interference in the experiment, a noise topology model is added, and parameter ] ′ (k) is augmented. erefore, (16) can be written as e j-th column in the model can be rewritten as SUAV Inputs SUAV outputs In a matrix form, (20) can be further rewritten as If the noise ] ′ (k) is a zero-mean unideal random series, the consistency and unbiased estimated value of θ j can be described as If j � 1, 2, . . . , m, the estimated values of the parameters in various rows, denoted as θ 1 , θ 2 , . . . , θ m , can be acquired. erefore, the estimated parameter value of the system can be acquired. e inverse calculation of a matrix demands heavy calculation, whereas a singular matrix easily appears during the test data processing. erefore, the following expressions can be set: e multivariable recursive augmented least squares identification algorithm process is shown in Figure 2 and can be consequently described as For the identification performance analysis, the following multi-input multioutput simulation model was established: where ](k) { } denotes the white noise N(0, 1) that follows a normal distribution in the diving phase. A four-displacement register with an input amplitude of 1 is used for the signal input.

Simulation Verification
To verify the effectiveness of the algorithm and analyze the identification performance of the modified multivariate augmented recursive least squares method, a simulation model is taken as an example. Figure 3 represents the comparison of the model noise identification results. e noise prediction of the system by the multivariate augmented recursive least squares method is closer to the real noise after approximation. e first channel noise prediction model characteristic parameters are within the range of ave n 1 , var n 1 and the second channel noise prediction model Mathematical Problems in Engineering 5 characteristic parameters are within the range of ave n 2 , var n 2 . e identification results are shown in Figure 4 and Table 1.

Flight Test
An SUAV generally follows a single-wing and V-shaped stern layout. Pitching and yawing are controlled by the V-shaped stern control surface while rolling is manipulated by the aileron. Electrically powered two-blade propellers are used for power supply, and the electric motor is installed behind the fuselage, i.e., a push-type SUAV is employed. As presented in Figure 5, a single ejection device is used for the launch. Table 2 lists the basic parameters of this SUAV. After launching the ejection device, the SUAV is controlled by the autopilot and a user in the loop to ensure that the SUAV flies in accordance with the preset trajectory at a    Mathematical Problems in Engineering   given velocity and height. During the flight, the target region can be monitored by an airborne camera. When a target is found, the SUAV dives, attacks, and destroys the target. Many measurement devices including GPS receiver, triaxial MEMS accelerometer, triaxial angular velocity gyroscope, air speedometer, and barometer are integrated into the navigation, guidance, and control system of SUAV. More specifically, longitude, latitude, height, and the heading information of the SUAV are provided by the GPS receiver. e accelerated velocity and angular velocity information are measured by the triaxial accelerometer and triaxial angular velocity gyroscope, respectively. e airspeed is measured by the airborne airspeed sensor.
Before the flight test, the UAV sensor and navigation equipment should be calibrated to reduce the noise impact. e flight test should be carried out in sunny weather and with low wind presence. Wind speed measurement curve of the flight test is displayed in Figure 6. e initial state of the drone was adjusted to the same settings before each dive test. To improve the validity of the data and reduce the effect of noise in the test data, data processing such as outlier processing, sampling time correction, data interpolation encryption, data coordination analysis, and data reconstruction were carried out on the preidentification data of the dive test [34,35]. e flight path of the SUAV is displayed in Figure 7. During testing, when the SUAV finds the target, it adjusts the attacking posture and consequently dives to attack. e SUAV performed five diving attacks in total. During the first four diving attacks, when the SUAV dived and approached the target, it was rapidly pulled up by the remote control handle to achieve attitude adjustment. During the fifth diving attack, the SUAV directly stroke the target. e adjusted pulse width outputs of each control manipulate channel during the simulated diving attack are displayed in Figure 8. e manipulated variables in the pitching and rolling channels exhibited the most drastic changes. Furthermore, these were consequently regarded as the main controlled variables that also constitute the coupled control input. e throttle channel during diving was controlled by the remote device under a constant input control pattern.

Results and Discussion
rough the model loss functions corresponding to different parameter combinations comparison, the model structure of the SUAV in the diving phase can be written as       where, A 1 , A 2 , B 0 , B 1 , B 2 , C 1 , and C 2 denote the parameter matrices to be identified. e first four sets of test data in the diving flight test were used as the identification data, while the fifth set of test data was used as the validation data. Table 3 lists the identification results. Figure 9 presents the convergence curves of the identification parameters of the matrices A 1 , A 2 , B 0 , B 1 , B 2 , C 1 , and C 2 . At the beginning of iteration, the identification parameters have significantly changed, which also increased the predicted model noise, as presented in Figure 10. Following approximately 100 iterations, the parameters exhibited gentle changes and slight adjustment. Figure 11 is a comparison curve of the identification results of a small drone dive model. e identification method proposed in this paper still predicts the changes in the pitch and roll angular rates well. e average errors of the pitch and roll angular rate are 0.0143 rad/s and 0.00868 rad/s, respectively, which is better than conventional least squares identification results (average error of 0.0251 rad/s and

Parameter
Identification value 0.00647 ave n ele 0.0101 var n ele 0.00841 ave n ail 0.00739 var n ail 0.0144 0.0125 rad/s). During the dive, the UAV's longitudinal motion range is larger, and the degree of nonlinearity is greater than that of the lateral direction, which makes the longitudinal identification accuracy slightly worse than that of the lateral direction.  conventional identification method (the average error is 0.1054 rad/s and 0.0902 rad/s).

Conclusions
Based on the flight test of the SUAV, the mathematical model of the dive phase was identified in this paper. A topological mathematical model of multi-input and multioutput noise terms for SUAV was constructed. Moreover, an improved multivariable augmented recursive least squares method was proposed. e above method was applied and analyzed via the identification test of a small UAV dive model. By comparing the results with the least squares method, it was shown that the proposed mathematical model and identification method can effectively obtain the parameters of the drone dive model. Furthermore, they can achieve high-precision identification of the drone dive model. is method is suitable for building nonmechanical modeling of the system when experimental data are available. In addition, the proposed method does not support the change mechanism of the internal physical parameters of the model. Future research work will focus on the application of tilt-rotor, trans-medium aircraft, and other new configuration aircraft tests. In addition, it will generate relevant test data, study the test modeling technology of special flight phases, and expand the scope of application of the proposed method. In terms of algorithm research, the results of more identification methods will be compared to continuously improve the identification performance of the model algorithm.

Data Availability
Some or all data, models, or codes generated or used during the study are proprietary or confidential in nature and may only be provided with restrictions (e.g., anonymized data).

Conflicts of Interest
e authors declare that they have no conflicts of interest.