## 4 Autonomous Landing Method
### Trajectory Generation and Tracking  
#### 4.3.1 Dynamic Model of the FPLIR  
In this section, two right-handed frames are used [48]--inertial frame \(\mathcal{F}  
The FPLIR is modeled as a multi-rotor robot, and then the dynamic model of position is expressed by:  
\[m\ddot{\boldsymbol{p}}^{\epsilon}= \begin{bmatrix}0\\ 0\\ -mg\end{bmatrix} +R_{b}^{\epsilon} \begin{bmatrix}0\\ 0\\ -T\end{bmatrix} +\boldsymbol{f}_{a} \tag{10}\]  
where \(m\) is the mass of the FPLIR, \(\boldsymbol{p}^{\epsilon}\) is the position of the FPLIR, \(g\) is gravitational acceleration, \(T\) is the total lift, and \(\boldsymbol{f}_{a}\) is the air damping force. The \(R_{b}^{\epsilon}\) is the rotation matrix from the north-east-down (NED) reference frame to the body frame (NED) and is given by:  
\[\boldsymbol{R}_{b}^{\epsilon}= \begin{bmatrix}\epsilon_{\theta}c_{\psi}&s_{ \theta}s_{\theta}c_{\psi}-c_{\phi}s_{\psi}&c_{\phi}s_{\theta}c_{\psi}+s_{\phi} s_{\psi}\\ c_{\theta}s_{\psi}&s_{\theta}s_{\theta}s_{\psi}+c_{\phi}c_{\psi}&c_{\phi}s_{ \theta}s_{\psi}-s_{\phi}c_{\psi}\\ -s_{\theta}&s_{\phi}c_{\theta}&c_{\phi}c_{\theta}\end{bmatrix}  \tag{11}\]  
with the shorthand notation abbreviated notations \(c_{x}\triangleq\cos x\), \(s_{x}\triangleq\sin x\), and \(t_{x}\triangleq\tan x\). The angular positions \(\varphi,\theta,\) and \(\psi\) denote, respectively, the roll, pitch, and yaw angle of the FPLIR. Therefore, Equation (11) can be written as:  
\[\boldsymbol{R}_{b}^{\epsilon}= \begin{bmatrix}c_{\theta}c_{\psi}&s_{\theta}s_{ \theta}c_{\phi}-c_{\theta}s_{\psi}&c_{\phi}s_{\theta}c_{\psi}+s_{\theta}s_{ \psi}\\ c_{\theta}s_{\psi}&s_{\theta}s_{\psi}+c_{\phi}c_{\psi}&c_{\phi}s_{\theta}s_{ \psi}-s_{\phi}c_{\psi}\\ -s_{\theta}&s_{\phi}c_{\phi}&c_{\phi}c_{\theta}\end{bmatrix}  \tag{12}\]  
\[ \begin{bmatrix}\ddot{p}_{n}\\ \ddot{p}_{e}\\ \ddot{p}_{d}\end{bmatrix} = \begin{bmatrix}0\\ 0\\ -g\end{bmatrix} +\frac{T}{m} \begin{bmatrix}c_{\phi}s_{\theta}c_{\psi}+s_{\theta} s_{\psi}\\ c_{\phi}s_{\theta}s_{\psi}-s_{\phi}c_{\psi}\\ c_{\phi}c_{\phi}\end{bmatrix} -\frac{k_{d}}{m} \begin{bmatrix}\dot{p}_{n}\\ \ddot{p}_{e}\\ \ddot{p}_{d}\end{bmatrix} \dot{p}_{e}\\ \ddot{p}_{d}\dot{p}_{d}\end{bmatrix} \tag{13}\]  
where \(p_{n}\), \(p_{e}\), and \(p_{d}\) represent positions in the north, east, and down directions, respectively. \(k_{d}\) is the friction coefficient. The air damping force is modeled as a proportional signed quadratic velocity. The dynamic equation of attitude is given by:  
\[\boldsymbol{I}^{b}\dot{\boldsymbol{\Omega}}^{b}+\boldsymbol{\Omega}^{b} \times\boldsymbol{I}^{b}\cdot\boldsymbol{\Omega}^{b}=\boldsymbol{\tau}+ \boldsymbol{g}_{a} \tag{14}\]  
where \(\boldsymbol{I}^{B}=diag\big{[}I_{xx}I_{yy}I_{zz}\big{]}\) is the rotational inertia, \(\boldsymbol{\Omega}^{b}= \begin{bmatrix}\Omega_{x}^{b}&\Omega_{y}^{b}&\Omega_ {y}^{b}\end{bmatrix} ^{T}\) is the vector of angular velocities in the body frame, \(\boldsymbol{\tau}= \begin{bmatrix}\tau_{x}&\tau_{y}&\tau_{z}\end{bmatrix} ^{T}\) is the control moment, and \(\boldsymbol{g}_{a}\) is the uncertainty of the model caused by higher aerodynamic effects and gravity shifts. Here, \(\boldsymbol{g}_{a}\) is given by:  
\[\boldsymbol{g}_{a}= \begin{bmatrix}s_{\alpha,\phi}\\ s_{\alpha,\theta}\\ s_{\alpha,\psi}\end{bmatrix} = \begin{bmatrix}J_{RP}\Omega_{y}^{b}(\varpi_{1}- \varpi_{2}+\varpi_{3}-\varpi_{4}+\varpi_{5}-\varpi_{6})\\ J_{RP}\Omega_{x}^{b}(-\varpi_{1}+\varpi_{2}-\varpi_{3}+\varpi_{4}-\varpi_{5}+ \varpi_{6})\\ 0\end{bmatrix}  \tag{15}\]  
Where \(J_{RP}\) is the rotational inertia of the motor rotor and propeller, and \(\varpi_{n}\) is the rotational speed of the motor. The components of angular velocity of the robot in the body frame are \(\Omega_{x}^{b}\), \(\Omega_{y}^{b}\), and \(\Omega_{z}^{b}\). These values are related to the derivatives of the roll, pitch, and yaw angles, which can be expressed as:  
\[ \begin{bmatrix}\dot{\boldsymbol{\phi}}\\ \dot{\theta}\\ \dot{\psi}\end{bmatrix} =