## 2 System and Problem Formulation
According to the above description, we have the following constraints:  
\[\Phi=\{\Theta[n],n=1,2,\ldots,N\} \tag{2}\] \[-\pi\leq\Theta[n]\leq\pi\] \[G[1]=(x_{g}[1],y_{g}[1],z_{g}[1])\] \[G[N]=(x_{g}[N],y_{g}[N],z_{g}[N])\] \[V_{min}\leq V[n]\leq V_{max},\ n=1,2,\ldots,N\] \[H_{l}\leq z_{g}[n]\leq H_{h}.\]  
The channel modeling of the system is measured by the channel coefficients at the \(n\)th time slot from the LS to PD via the RIS as \(h[n]\), which are represented by (Note: Here, due to the turbulence of atmosphere, is a random variable that only relies on the constant from  
Figure 1: The position and orientation of the system in the considered coordinates. Figure 2: The visual LS and the project center of laser. circumstances, most of the papers engage the UAV and FSO communication optimization and do not take it into consideration, such as [38; 39])  
\[\Theta[n](\theta_{1}[n],\theta_{2}[n])\in(-\pi,\pi). \tag{1}\]0  
where \(h_{t}\) is the turbulence of the atmosphere, this parameter is related to the wind speed in refractive index and the laser wavelength [40], where the probability density function (PDF) of \(h_{t}\) can be denoted as [41]  
\[\Theta[n](\theta_{1}[n],\theta_{2}[n])\in(-\pi,\pi). \tag{1}\]1  
where \(K_{\alpha-\beta}(. )\) is the modified Bessel function of the second kind, the \(\Gamma(. )\) is the Gamma function, and the \(\alpha\) and \(\beta\) denote the numbers of small and large turbulence cells as follows,  
\[\Theta[n](\theta_{1}[n],\theta_{2}[n])\in(-\pi,\pi). \tag{1}\]2 \[\Theta[n](\theta_{1}[n],\theta_{2}[n])\in(-\pi,\pi). \tag{1}\]3  
respectively, where Rytov variance \(\sigma_{R}\), which can be calculated as \(\sigma_{R}^{2}=0.5k^{\frac{\alpha}{2}}C_{N}^{2}L^{\frac{11}{2}}\), and \(C_{N}^{2}\) is the index of refraction structure parameter, \(k=\frac{2\pi}{\lambda_{f}}\), can be obtained by optical wavelength \(\lambda_{f}\). The \(h_{a}[n]\) is the atmospheric loss caused by fog, which can be given by  
\[\Theta[n](\theta_{1}[n],\theta_{2}[n])\in(-\pi,\pi). \tag{1}\]4  
where \(\sigma_{c}\) and \(\sigma_{f}\) denote the measurement of atmospheric attenuation conditions under clean air and fog [42], respectively, which can be given as \(\sigma_{c}=\frac{3.91}{L_{v}}\left(\frac{\lambda}{550}\right)^{-p_{c}}\) and \(\sigma_{f}=\frac{3.91}{L_{v}}\left(\frac{\lambda}{550}\right)^{-p_{f}}\) where  
\[p_{c}=\left\{ \begin{array}{cc}1.6&L_{v}>50\\ 1.3&6<L_{v}<50,\end{array} \right. \tag{7}\]  
and  
\[p_{f}=\left\{ \begin{array}{cc}0.16L_{v}+0.34&1<L_{v}<6\\ L_{v}-0.5&0.5<L_{v}<1\\ 0&L_{v}<0.5,\end{array} \right. \tag{8}\]  
where \(L_{v}\) denotes the visibility range in air. Furthermore, \(d_{c}[n]\) and \(d_{f}[n]\) denote the laser propagation distances in clear air and fog at the \(n\)th time slot, respectively, which can be obtained as in Appendix A. Due to the UAV operating with a 3D flight, the incident laser beam is not perpendicular to the PD, and the Gaussian profile beam footprint is an ellipse on the PD plane as shown in Figure 3. The pointing error loss \(h_{p}[n]\) is caused by laser power fraction on the PD and can be the integral results from the beam power intensity, denoted as  
\[I(\mathbf{\rho};[n])=\frac{2}{\pi W_{x}[n]W_{y}[n]}\exp\left(-\frac{2||\mathbf{\rho}||^ {2}}{W_{x}[n]W_{y}[n]}\right), \tag{9}\]then consider it under the polar coordinate system of the detector plane with radial \(a\) and angular \(2\pi\), we can express the pointing error loss as  
\[ \begin{split} h_{p}[n]=&\int_{A}I(||\mathbf{\rho}[n]-\mathbf{r} [n]||)\mathrm{d}\mathbf{\rho}[n]\\ =&\int_{0}^{a}\int_{0}^{2\pi}\frac{2}{\pi W_{y}[n]W_{ x}[n]}\exp\frac{-2(\mathbf{\rho}[n]-\mathbf{r}[n])^{2}}{W_{y}[n]W_{x}[n]}\mathbf{\rho}[n] \mathrm{d}\mathbf{\phi}[n]\mathrm{d}\mathbf{\rho}[n]\end{split}