## 2 Skywave Ionosphere Channel for Hypersonic Vehicles
\begin{array}{l}\frac{\partial X}{\partial r}\approx\frac{80.8}{f^{2}}\frac{ N_{e}(r+\Delta r,\theta,\varphi)-N_{e}(r,\theta,\varphi)}{\Delta r}\\ \frac{\partial X}{\partial\theta}\approx\frac{80.8}{f^{2}}\frac{N_{e}(r, \theta+\Delta\theta,\varphi)-N_{e}(r,\theta,\varphi)}{\Delta\varphi}\\ \frac{\partial X}{\partial\varphi}\approx\frac{80.8}{f^{2}}\frac{N_{e}(r, \theta,\varphi+\varphi)-N_{e}(r,\theta,\varphi)}{\Delta\varphi}\end{array}  \tag{9}\]  
Based on the initial \(r\), \(\theta\), \(\varphi\), and \(\alpha\), the ray equation can be solved iteratively via the step size of the group path \(P^{\prime}\). This allows for the calculation of \(r\), \(\theta\), \(\varphi\), \(k_{r}\), \(k_{\varphi}\), and \(k_{\theta}\) for each point along the ray path, and, by connecting these points, the trajectory of the skywave propagation can be obtained in a spherical coordinate system. If electromagnetic waves are emitted from the ground, the latitude and longitude of the emission point are denoted by \(\lambda_{0}\) and \(\varphi_{0}\), respectively. The launch angle starts from the south direction of the ground, with the positive direction pointing towards the zenith, denoted as \(\beta_{0}\). The launch azimuth starts from the northern direction on the ground and is considered positive when measured clockwise, denoted as \(\alpha_{0}\). The initial values for the Haselgrove equation system for rays are set as follows:  
\[\left\{ \begin{array}{l}r_{0}=R_{0}\\ \theta_{0}=\frac{\pi}{2}-\lambda_{0}\\ \varphi_{0}=\varphi_{0}\\ \textbf{{k}}_{r0}=\textbf{{k}}\cos(\frac{\pi}{2}-\beta)=\frac{\omega}{c}\sin \beta_{0}\\ \textbf{{k}}_{\theta 0}=-\textbf{{k}}\cos\beta\cos\alpha=-\frac{\omega}{c}\cos \beta_{0}\cos\alpha_{0}\\ \textbf{{k}}_{\varphi 0}=\textbf{{k}}\cos\beta\cos(\frac{\pi}{2}-\alpha)=\frac{ \omega}{c}\cos\beta_{0}\sin\alpha_{0}\end{array} \right. \tag{10}\]  
where \(R_{0}\) is the radius of the Earth, typically taken as 6371 km. As shown in Figure 2, let (\(r_{M}\), \(\theta_{M}\), \(\varphi_{M}\)) be the \(k\)-th ray that reaches the current position of the vehicle, where \(M\) represents the total number of iterations for group path \(P^{\prime}\) at time \(t\). In Formula (6), \(\int\alpha_{iok}(t)dl(t)\) can be written as  
\[\int\alpha_{iok}(t)dl(t)=\sum_{1}^{i=M}\omega\cdot sqrt\left(\frac{\mu_{0} \varepsilon_{0}}{2}\left(\sqrt{\varepsilon_{r}^{2}(r_{i},\theta_{i},\varphi_{ i})+\left(60\lambda\sigma(r_{i},\theta_{i},\varphi_{i})\right)^{2}}-\varepsilon_{r}(r_{i}, \theta_{i},\varphi_{i})\right)\right)\cdot d^{\prime}(r_{i},\theta_{i},\varphi_ {i}) \tag{11}\]  
where \(l\) (\(t\)) is the ray path traveled by the electromagnetic wave from the emission point to the vehicle. The conductivity and equivalent relative dielectric constant of the ionosphere at this location are \(\sigma\) (\(r_{i}\), \(\theta_{i}\), \(\varphi_{i}\)) and \(\varepsilon_{r}\) (\(r_{i}\), \(\theta_{i}\), \(\varphi_{i}\)). \(d^{\prime}\) (\(r_{i}\), \(\theta_{i}\), \(\varphi_{i}\)) is the distance between the previous ray trajectory position point (\(r_{i-1}\), \(\theta_{i-1}\), \(\varphi_{i-1}\)) and this ray trajectory position point (\(r_{i}\), \(\theta_{i}\), \(\varphi_{i}\)):  
\[d^{\prime}(r_{i},\theta_{i},\varphi_{i})=\sqrt{r_{i}^{2}+r_{i-1}^{2}-2r_{i}r_{ i-1}[(\sin\theta_{i}\sin\theta_{i-1}\cos(\varphi_{i}-\varphi_{i-1})+\cos \theta_{i}\cos\theta_{i-1})]} \tag{12}\]  
For skywave OTH channels, the Haselgrove ray equation is used for iterative propagation within the medium. If the ray path reaches the ground, the iteration terminates, and the complete ray path cannot be fully estimated. Typically, as the distance between the vehicle and the launch station increases, the ray (signal) experiences multiple reflections between the ionosphere and the ground before reaching the vehicle. Therefore, updating the elevation angle \(\beta_{u}\) of the ground-reflected ray during the iteration process is crucial, where \(u\) is the number of times the ray is reflected to the ground.