# Ionospheric scintillation modeling II. Theoretical study of power spectral densities of ionospheric scintillation

Grégory Morel, Vincent Fabbro, Olivier Boisot, Laurent Feral

► To cite this version:

Grégory Morel, Vincent Fabbro, Olivier Boisot, Laurent Feral. Ionospheric scintillation modeling II. Theoretical study of power spectral densities of ionospheric scintillation. Journal of Space Weather and Space Climate, 2025, 15, pp.54. 10.1051/swsc/2025050. hal-05444639

HAL Id: hal-05444639

https://hal.science/hal-05444639v1

Submitted on 6 Jan 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Ionospheric scintillation modeling II. Theoretical study of power spectral densities of ionospheric scintillation

Grégory Morel<fnref n="1" />*, Vincent Fabbro<fnref n="1" />, Olivier Boisot<fnref n="2" />, and Laurent Féral<fnref n="1" />

Received 23 April 2025 / Accepted 16 October 2025

**Abstract**— An electromagnetic signal that has passed through a layer of ionospheric irregularities can undergo significant modifications. The Power Spectral Density (PSD) of such a signal can then be rich in information about the medium through which it has passed. Thus, this paper is dedicated to the theoretical study of the PSD of the log-amplitude and phase of a signal received by a ground receiver, for a satellite-Earth link that encounters ionospheric irregularities. Incident plane wave (PW), incident spherical wave (SW) and corrected plane wave (CPW) formalisms are considered. For each approach, an analytic integral formulation of the log-amplitude and phase PSD of the received signal are established and compared under the assumptions of weak scattering and a thin layer of irregularities relative to the length of the satellite-Earth link. The SW analytic formalism is validated trough comparison with a numerical approach solving the Parabolic Wave Equation (PWE) in spherical coordinates combined with Multiple Curved Phase Screens (MCPS). Asymptotic derivations of the PSD analytic formulations for low and high frequencies are proposed, and in addition, particular attention is given to the frequency corresponding to the maximum of the log-amplitude PSD. These formulations are analysed, and their sensitivity to the characteristics of the medium is discussed to address the inversion problem, consisting of deriving parameters that describe the irregularities of the medium through its log-amplitude and phase PSD.

**Keywords:** Propagation / Ionospheric scintillation / Taylor hypothesis / Parabolic Wave Equation / Power Spectral Density / Fresnel scale / Rytov theory / Weak scattering

## 1 Introduction

This paper is the second part of a study on ionospheric scintillation modeling and analysis. It addresses the power spectral analysis (PSD) of ionospheric scintillation in a plane wave (PW) and a spherical wave (SW) approach. The first part focuses on the analysis of scintillation indices and is the topic of a companion paper.

Random variability of the received electromagnetic field due to fluctuations of electronic density in the ionosphere can cause degradation in positioning and communication system performances. Thus, accurate knowledge of the amplitude and phase signal received can help data processing research to mitigate these effects. They are related to the diffusion of electromagnetic waves in the ionosphere, which can present strong turbulent irregularities. Characterizing these ionospheric irregularities and their impact on the propagated field remain a challenge. However, the theoretical interpretation of the variability of the

received signal can provide information about the medium through which the waves have propagated.

Nowadays, many two-frequency radionavigation systems allow for Total Electron Content (TEC), amplitude and phase scintillation measurements of the signal. The set of these measurements potentially contains a wide range of information describing the ionospheric medium and its variability. Mastering an analytic formalism to describe signal variations in response to medium fluctuations is essential for tackling an inversion problem effectively. In addition to understanding the current ionospheric conditions, mastering the long-term climatology of amplitude and phase scintillation phenomena allows a better design of electromagnetic systems by associating an operational availability. Statistics of scintillation events are therefore relevant. These are linked to the long-term solar activity and its cycle, the Earth's magnetic activity, the season, the local time, the position on the Earth and the geometry of the radio-electric link with respect to the Earth's magnetic field.

In the past, combined observations of incoherent scatter and scintillation observations have led to a better understanding of the morphology of ionospheric irregularities (Aarons, 1982;

<footnote n="1">¹ DEMR, ONERA, Université de Toulouse, 31000 Toulouse, France</footnote>

<footnote n="2">² DEMR, ONERA, 13330 Salon-de-Provence, France</footnote>

<footnote>*Corresponding author: gregory.morel@onera.fr</footnote>
Livingston et al., 1982; Rino et al., 1983). In this context, the theoretical work of modelling the scintillation indices in amplitude, S₄ and phase, σφ, were particularly useful (Rino, 1979). More recently, in Carrano et al. (2016) the authors propose a method to estimate the zonal irregularity drift from measurements of amplitude and phase scintillation indices by a single GNSS station. The Fresnel frequency fF has also been studied by several authors to characterize ionospheric medium irregularities (McCaffrey & Jayachandran, 2019; Song et al., 2022; Sun et al., 2023). It is given by fF = V/(√2rF) with V the relative layer drift (composed by ionospheric piercing point of the LOS and the irregularity drift), rF = √λLv the Fresnel distance, where λ is the wavelength and Lv the ionospheric layer-receiver distance. In the literature, the Fresnel frequency is sometimes associated with the frequency of the maximum of the log-amplitude PSD (that is also called rollover frequency) (Forte & Radicella, 2002). In this paper, the limits of this assumption are specifically studied. Going forward, to describe morphological and spectral features of irregularities, inversion of GPS observations has been used. For instance, in Conroy et al. (2021), Vaggu et al. (2023) an inverse method is applied to derive the parameters that describe the ionospheric layer irregularities. The approach estimates the best fit of numerical model outputs to GPS observations. The numerical model is a 3D-PWE/2D-MCPS (Parabolic Wave Equation method with Multiple Curved Phase Screens) resolution of the Helmholtz propagation equation, and the inversion approach is based on spectral analysis. This approach requires the use of a numerical regression algorithm. To optimize computation times and facilitate the inversion of results, the development of an analytic formalism is fundamental.

The modelling of the amplitude and phase of the received signal has been the subject of intensive theoretical research and several approaches can be found in the open literature (Yeh & Liu, 1982; Vasylyev et al., 2022). Notably, one can cite the work of Tatarskii (1971), Rytov et al. (1989) and Wheelon (2003), which led to an asymptotic approach under low disturbance and thin layer hypothesis. One strong feature of these formulations is that they would allow for attempts to access certain parameters, describing the statistics of variability of the ionospheric environment.

The originality of this paper, by comparison to Fabbro & Feral (2012) and Galiègue (2015), is to develop analytic formulations of log-amplitude and phase PSD, which could be used to make inversion of ionospheric parameters describing irregularities. To facilitate the inversion and interpretation of spectra, asymptotic formulations at low and high frequencies are also proposed, as well as an approximate formulation of the maximum of the log-amplitude spectrum. The goal is to provide to the community key formulations or parameters, giving more direct access to the description of the ionospheric turbulence. The analysis is carried out in terms of log-amplitude and phase Power Spectral Densities (PSD). To validate the presented results, both 3D-PWE/2D-MCPS numerical scheme and Rytov's theory of smooth perturbation results for spherical waves are compared, and differences between PSD in plane and spherical wave formalisms are also studied. In addition, analytic derivations of the frequency of the maximum of the log-amplitude PSD obtained under the weak scattering and thin-layer hypotheses are also proposed.

Knowing the spectral characteristics of the ionospheric medium makes it possible to better describe the dynamics of the scintillation of the received signal, to describe its anisotropy, to understand its cascade of energy in the wavenumber domain by studying its inertial slope, and to know the size of irregularities via its outer scale of turbulence. As the measured scintillation, through the amplitude and phase scintillation indices (S₄ and σφ) studied in the first paper, have different dynamics of variation and different statistical occurrence depending on the geophysical conditions, we can then imagine different characteristics of the PSD of the irregularities depending on these conditions. Knowledge of this local climatology can be particularly relevant for system design.

The paper is organized as follows: Section 2 is focused on analytic formulations of the ionospheric scintillation impact in terms of log-amplitude and phase PSD for both incident Plane Wave (PW) and Spherical Wave (SW). In addition to these formulations, a heuristic approach, called Corrected Plane Wave (CPW) is proposed as a modification of the plane wave formalism. The frequency of the maximum of the log-amplitude spectrum is then derived, and its accuracy is studied. In Section 3, the different approaches are compared (i.e., analytic incident PW/SW/CPW formulations), additionally, the sensitivity of the spectra to ionospheric parameters is analyzed. Finally, the different results obtained are summarized and discussed in Section 4 which concludes the paper.

## 2 PSD analytic formulations under weak scattering assumption

In Galiègue et al. (2016), the spectrum of the electronic density fluctuations SΔNₐ has been derived in the line of sight coordinate system (u, v, s). It is given by:

$$S_{\Delta N_e}(k_u, k_v, k_s) = a_x^{3-p_m} A_y A_z C_s$$

$$\left( A k_u^2 + B k_v^2 + 2 C k_u k_v + f(k_s) + \frac{K_0^2}{a_x^2} \right)^{-p_m/2}, \quad (1)$$

where A, B, and C are the anisotropic coefficients expressed in the line of sight basis (cf. expressions in Galiègue et al., 2016). They depend on the link geometry: γ, αz and ψ angles (see Fig. 1 of Morel et al., 2025) and anisotropic ratios A_y, A_z, a_x is a dimensionless scaling factor, K₀ = 2π/L₀ [rad m⁻¹] is the wave number associated with the outer scale of the turbulence and p_m the slope of the Power Spectral Density of the irregularities, finally C_s [m⁻³-p_m] is the strength of the ionospheric turbulence (Rino, 2011; Galiègue et al., 2016).

From this description, the Rytov's theory of smooth perturbations (Rytov et al., 1989; Wheelon, 2003), allows to derive analytic formulations of log-amplitude and phase spectra. The electric field E₁(R) can then be decomposed in terms of amplitude A(R) and phase φ(R) with E₁(R) = A(R)e^(i(φ₀+φ(R))). For more details, the reader should refer to Morel et al. (2025).

The signal first traverses a slant distance L_t from the transmitter to the upper part of the ionospheric layer, then covers a distance Rᵢonₐ within the ionospheric layer, and finally travels a distance L_v from the lower part of the ionospheric layer to the ground station. The overall distance traveled by the signal
is $R = L_t + R_{\text{iono}} + L_v$. In this context, H corresponds to the altitude of the ionospheric layer, and $\Delta H$ indicates its thickness as described in more detail in Morel et al. (2025).

## 2.1 Log-amplitude and phase PSD

The log-amplitude $W_\chi$ and phase $W_\varphi$ spectra are respectively the temporal Fourier transform of the log-amplitude and phase covariances $\langle\chi(\mathbf{R}, t)\chi(\mathbf{R}, t+\tau)\rangle$, and $\langle\varphi(\mathbf{R}, t)\varphi(\mathbf{R}, t+\tau)\rangle$, assuming that both processes are stationary (Fabbro & Feral, 2012):

$$W_\chi(\omega) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} \langle \chi(\mathbf{R}, t) \chi(\mathbf{R}, t+\tau) \rangle e^{-j\omega\tau} d\tau, \quad (2)$$

$$W_\varphi(\omega) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} \langle \varphi(\mathbf{R}, t) \varphi(\mathbf{R}, t+\tau) \rangle e^{-j\omega\tau} d\tau, \quad (3)$$

where $\omega$ [rad s⁻¹] is the pulsation of the propagated wave, $t$ and $\tau$ [s] re the time, ($\tau$ is the dual variable of $\omega$ in the Fourier space). The transformation between space and time is introduced by assuming Taylor's hypothesis (Taylor, 1938): the medium is considered frozen and advected with a uniform speed **V** (Yeh & Liu, 1982). This velocity vector **V** is defined by its coordinates in the Line Of Sight (LOS) basis **V** = ($V_u, V_v, V_s$), which allows defining the drift speed transverse to the line of sight $V_r = \sqrt{V_u^2 + V_v^2}$. In such conditions, it can be shown that (Ishimaru, 1972; Wheelon, 2003):

$$\langle \chi(\mathbf{R}, t) \chi(\mathbf{R}, t+\tau) \rangle = \pi r_e^2 \lambda^2 R_{\text{iono}} \int_{-\infty}^{+\infty} \int_{-\infty}^{+\infty} S_{\Delta N_e}(k_u, k_v, 0) \mathcal{F}_\chi(k_u, k_v) e^{\int (k_u V_u + k_v V_v) dt} dk_u dk_v, \quad (4)$$

$$\langle \varphi(\mathbf{R}, t) \varphi(\mathbf{R}, t+\tau) \rangle = \pi r_e^2 \lambda^2 R_{\text{iono}} \int_{-\infty}^{+\infty} \int_{-\infty}^{+\infty} S_{\Delta N_e}(k_u, k_v, 0) \mathcal{F}_\varphi(k_u, k_v) e^{\int (k_u V_u + k_v V_v) dt} dk_u dk_v, \quad (5)$$

where $F_\chi$ and $F_\varphi$ are the filtering functions in amplitude and phase respectively (Morel et al., 2025). By combining equations (2) and (4), and equations (3) and (5), it comes (Ishimaru, 1972; Wheelon, 2003; Fabbro & Feral, 2012):

$$W_\chi(\omega) = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}}}{V_u} \int_{-\infty}^{+\infty} S_{\Delta N_e} \left( \frac{\omega - k_v V_v}{V_u}, k_v, 0 \right) \mathcal{F}_\chi \left( \frac{\omega - k_v V_v}{V_u}, k_v \right) dk_v, \quad (6)$$

$$W_\varphi(\omega) = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}}}{V_u} \int_{-\infty}^{+\infty} S_{\Delta N_e} \left( \frac{\omega - k_v V_v}{V_u}, k_v, 0 \right) \mathcal{F}_\varphi \left( \frac{\omega - k_v V_v}{V_u}, k_v \right) dk_v. \quad (7)$$

In these formulations, the incident wave illuminating the ionosphere irregularites can be a spherical wave (SW) or a plane plane (PW). SW or PW specificity can be taken into account by changing the spectral filtering function, named $\mathcal{F}_{\chi,\varphi}^{\text{PW}}$ and

$\mathcal{F}_{\chi,\varphi}^{\text{SW}}$, as realized in (Morel et al., 2025) to compute the amplitude and phase variances. The plane wave formulation of the filtering functions are:

$$\begin{aligned} \mathcal{F}_{\chi, \varphi}^{\text{PW}}(k_u, k_v) = & 1 \mp \text{sinc}\left(\frac{(k_u^2 + k_v^2) r_F^{\text{PW}2} R_{\text{iono}}}{4\pi L_v}\right) \\ & \cos\left(\frac{(k_u^2 + k_v^2) r_F^{\text{PW}2}}{2\pi}\left(1 + \frac{R_{\text{iono}}}{2L_v}\right)\right), \end{aligned} \quad (8)$$

with $r_F^{\text{PW}} = \sqrt{\lambda L_v}$ the Fresnel radius for the incident plane wave. For the spherical wave formulation, (8) becomes:

$$\begin{aligned} \mathcal{F}_{\chi, \varphi}^{\text{SW}}(k_u, k_v) = & 1 \mp \frac{2}{R_{\text{iono}}} \sqrt{\frac{k_0(L_t + R_{\text{iono}} + L_v)}{4(k_u^2 + k_v^2)}} \times \\ & \left[ \cos\left(\frac{(L_t + R_{\text{iono}} + L_v)(k_u^2 + k_v^2)}{4k_0}\right)(C(t_2) - C(t_1)), \right. \\ & \left. + \sin\left(\frac{(L_t + R_{\text{iono}} + L_v)(k_u^2 + k_v^2)}{4k_0}\right)(S(t_2) - S(t_1)) \right] \end{aligned} \quad (9)$$

where $t_1$ and $t_2$ are given in the first paper ans C() and S() are the Fresnel integrals. This way, similarly to variance formulation correction in Morel et al. (2025), it is here proposed to correct the plane wave filtering functions by substitution of the Fresnel distance to get $\mathcal{F}_{\chi,\varphi}^{\text{CPW}}$. The filtering functions are then given by:

$$\begin{aligned} \mathcal{F}_{\chi, \varphi}^{\text{CPW}}(k_u, k_v) = & 1 \mp \text{sinc}\left(\frac{(k_u^2 + k_v^2) r_F^{\text{SW}2} R_{\text{iono}}}{4\pi L_v}\right) \\ & \cos\left(\frac{(k_u^2 + k_v^2) r_F^{\text{SW}2}}{2\pi}\left(1 + \frac{R_{\text{iono}}}{2L_v}\right)\right), \end{aligned} \quad (10)$$

with $r_F^{\text{SW}} = \sqrt{\frac{\lambda L_v}{L_v + L_t}}$, the Fresnel radius for an incident spherical wave. The impact of the incident plane wave is analyzed in more details in Section 3.

## 2.2 Log-amplitude and phase PSD asymptotes

Whatever the filtering function (i.e., under SW, PW, or CPW approximation), the analytic expressions of the integrals in equations (6) and (7), are complicated to derive. However, it is possible to obtain asymptotic expressions at low and high frequencies that can be useful for inversion purposes. Assuming a thin turbulent layer compared to the layer-receiver distance ($R_{\text{iono}} \ll L_v$) and considering simplified weighting functions, asymptotic formulations can be derived for the log-amplitude and phase spectra. The details of the mathematical derivations are given in Appendix A. They lead to the following high-frequency behavior of log-amplitude and phase PSD (similar for plane wave and spherical wave formalisms):

$$\begin{aligned} W_{\chi, \varphi}^{\text{PW,SW}}(\omega \to +\infty) = & \frac{\pi \sqrt{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s \Gamma(\frac{p_m-1}{2})}}{\Gamma(p_m/2)(AV_v^2 + BV_u^2 - 2CV_u V_v)^{1-p_m/2}} \times \\ & \left(\omega^2(AB - C^2) + \frac{K_0^2}{a_x^2}(AV_v^2 + BV_u^2 - 2CV_u V_v)\right)^{\frac{1-p_m}{2}}. \end{aligned} \quad (11)$$
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Comparison of theoretical and numerical PSD for SW formalism. (a) Log-amplitude PSD for SW formalism. (b) Phase PSD for SW formalism.</figcaption></figure>

## 2.4 Maximum of the log-amplitude PSD

Different authors consider that the peak value of the log-amplitude spectrum, corresponding to its maximum value as a function of frequency, occurs for the Fresnel frequency $f_F$ (Yeh & Liu, 1982; Forte & Radicella, 2002; Wernik et al., 2003; McCaffrey & Jayachandran, 2019). However, this is a hypothesis that deserves to be quantified. To do this, an exercise of estimating the difference between the frequency of the PSD maximum and the Fresnel frequency was performed. Note that the Fresnel frequency is approximately the first maximum of the filtering function for log-amplitude. Indeed, for the analysis, it is interesting to compute the value of $k_r = \sqrt{k_u^2 + k_v^2}$ corresponding to the first null of the cardinal sine and cosine:

$$k_{\sin c} = \frac{2\pi}{\sqrt{\lambda R_{\text{iono}}}}, \qquad (15)$$

$$k_{\cos} = \frac{2\pi}{\sqrt{4\lambda L_v \left(1 + \frac{R_{\text{iono}}}{2L_v}\right)}}. \qquad (16)$$

With the thin layer assumption $R_{\text{iono}} \ll L_v$, it follows that $k_{\cos} \ll k_{\sin c}$, so the first maximum value of the filtering functions corresponds approximately to the first $-1$ value of the cosine term, i.e., for:

$$k_m \approx \frac{2\pi}{\sqrt{2} r_F^{\text{PW,SW}}}. \qquad (17)$$

Finally, with the relationship between frequency and wave number, the Fresnel frequency can be expressed as:

$$f_F^{\text{PW,SW}} = \frac{V_r}{\sqrt{2} r_F^{\text{PW,SW}}}. \qquad (18)$$

For the simulations, two different configurations have been defined, corresponding to a LEO satellite moving over the polar

or equatorial region, where ionospheric scintillation phenomena are frequent. The parameters describing these configurations are reported in Table 1. The main difference between polar and equatorial conditions is the orientation of the Earth's magnetic field parametrized by angles $(\gamma, \psi)$ and the anisotropy of the ionospheric irregularities, characterized by the anisotropic ratios $(A_y, A_z)$. The drift velocity of the medium is set to $1000 \text{ m s}^{-1}$ for the polar configuration and $100 \text{ m s}^{-1}$ for the equatorial configuration. However, the velocity of the ionospheric piercing point must also be taken into account, as explained in (Morel et al., 2025). Naturally, this velocity depends on the satellite's altitude. The other parameters have been kept equal between the two configurations, to facilitate the comparison. Note that the satellite altitude is only used for the SW formalism.

Figure 2 shows the relative error between the frequency of the maximum of the temporal log-amplitude spectrum (eq. (6)) estimated numerically and the Fresnel frequency (eq. (18)) for the polar (Fig. 2a) and equatorial (Fig. 2b) configurations. The study is proposed according to the main anisotropy $A_z$ and the slope variations of PSD $p_m$, corresponding to the abscissa and ordinate axes, respectively, whereas the relative error is reported in a color scale. Because relative error is plotted, results are the same regardless of the approach used, PW or SW.

Only for values of slope lower than 4 and high values of anisotropy, the Fresnel frequency seems to roughly approximate the frequency of the maximum $f_{\text{max}}$, though the relative error is still greater than 10% even in the most favorable cases.

So the statement that the Fresnel frequency corresponds to the frequency of the maximum of the log-amplitude spectrum is an approximation and does not really correspond to the frequency of the maximum of the log-amplitude spectrum in all configurations, as shown in Figure 2. A better expression for the frequency of the maximum has been researched here. Initially, we aim to find $\omega_{\text{max}}$ such that $\left[ \frac{dW_{\gamma}(\omega)}{d\omega} \right]_{\omega=\omega_{\text{max}}} = 0$, but since expression (6) makes an integral appear, obtaining
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Relative error in percentage between the frequency of the maximum of the log-amplitude spectrum estimated numerically from (6) and the Fresnel frequency (18), $|f_{\max}-f_F|$, in color scale for both polar (a) and equatorial (b) configurations according to the anisotropic ratio $A_z$ and the slope of the irregularities PSD $p_m$.</figcaption></figure>

an analytic expression for $\omega_{\max}$ is extremely challenging. The key idea is therefore to first determine $\mathbf{k}_{\max}$ such that $\nabla W_{\chi}(\mathbf{k}) = (0,0)$, where $W_{\chi}(\mathbf{k}) = W_{\chi}(k_u, k_v)$ is the 2D spatial spectrum given by:

$$W_{\chi}(k_u, k_v) = \pi r_e^2 \lambda^2 R_{\text{iono}} S_{\Delta N_e}(k_u, k_v, 0) F_{\chi}(k_u, k_v), \quad (19)$$

where the demonstration for this derivation can be found in Wheelon (2003), Morel et al. (2024) for an incident plane wave and an incident spherical wave, respectively. The expression of the 2D spatial spectrum does not make an integral appear as the temporal spectrum defined by equation (6). It is therefore easier to calculate the spectrum derivative and its roots to estimate the maximum. The mathematical details are reported in Appendix B.1. Two methods are used, providing two different solutions:

$$k_{\max} = \begin{cases} \sqrt{\frac{4K_0^2}{p_m - 4}} & \text{if } p_m > 4 + 2.3\Delta \\ \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)}} & \text{otherwise,} \end{cases} \quad (20)$$

with $\Delta = \frac{K_0^2 r_F^2}{2\pi}$, $a_x$ and $G$ do not appear here because we are in the particular case where $A_y = A_x = 1$ and so $a_x = G = 1$. The transformation from the spatial domain to the temporal domain is achieved using Taylor's hypothesis, allowing wave numbers to be converted into frequencies. It is possible to find a relation between the temporal spectrum and the spatial one by using the following relation $W_{\chi}(\mathbf{k}, \omega) = \delta(\mathbf{k} \cdot \mathbf{V}_r - \omega)W_{\chi}(\mathbf{k})$, which can be derived from (Tatarskii, 1971).

$$W_{\chi}(\omega) = \int_{-\infty}^{+\infty} W(\mathbf{k}, \omega) d^2\mathbf{k}, \quad (21)$$

$$= \int_{-\infty}^{+\infty} \delta(\mathbf{k} \cdot \mathbf{V}_r - \omega) W_{\chi}(\mathbf{k}) d^2\mathbf{k}, \quad (22)$$

$$= \int_{-\infty}^{+\infty} \delta(k_u V_u + k_v V_v - \omega) W_{\chi}(k_u, k_v) d^2\mathbf{k}, \quad (23)$$

$$= \frac{1}{V_u} \int_{-\infty}^{+\infty} \delta\left(k_u - \frac{\omega - k_v V_v}{V_u}\right) W_{\chi}(k_u, k_v) d^2\mathbf{k}, \quad (24)$$

$$= \frac{1}{V_u} \int_{-\infty}^{+\infty} W_{\chi}\left(\frac{\omega - k_v V_v}{V_u}, k_v\right) dk_v. \quad (25)$$

This last relationship can be verified by independently computing temporal (Eq. (6)) and spatial (Eq. (19)) spectra leading to equality (25). Let $\theta$ be the angle between the u-axis and the principal direction of the 2D spectrum, i.e., where the energy is mainly concentrated, $k_{\max}$ can be decomposed as: $\mathbf{k}_{\max} = k_{\max}\cos(\theta)\mathbf{u} + k_{\max}\sin(\theta)\mathbf{v}$. Using now the dispersion relationship $\omega = 2\pi f = k_u V_u + k_v V_v$, the frequency $f_{\max}$ corresponding to the wave number $k_{\max}$ can be derived as:

$$f_{\max} = \left| \frac{1}{2\pi} (k_{\max} \cos(\theta) V_u + k_{\max} \sin(\theta) V_v) \right|. \quad (26)$$

However, it is important to keep in mind that this transformation is an approximation, as the frequency corresponding to the maximum of the temporal log-amplitude spectrum cannot be determined analytically due to the integral in equation (6). In particular cases where the Earth's magnetic field $\mathbf{H}_0$ (z axis) is aligned with u or v, or if one of the two velocity components is zero, the frequency of the maximum is:

$$f_{\max} = \begin{cases} \frac{V_r}{2\pi} \sqrt{\frac{4K_0^2}{p_m - 4}} & \text{if } p_m > 4 + 2.3\Delta \quad (27a) \\ \frac{V_r}{2\pi} \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)}} & \text{otherwise,} \end{cases}$$

with $V_r = V_u$ or $V_v$. An attempt to improve the accuracy was made by applying one iteration of the Newton-Raphson
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Relative error in percentage between the frequency of the maximum of the log-amplitude spectrum estimated numerically from (6) and the estimated one from equations (27a) and (27b) in color scale for both polar (a) and equatorial (b) configurations according to the anisotropic ratio $A_z$ and the slope of the irregularities PSD $p_m$.</figcaption></figure>

algorithm. The gain is not significant, and the obtained analytic expressions are more complicated, and thus would be difficult to use for inversion parameters. The results of this analysis are, however, reported in Appendix B.1, for information.

## 2.5 Assessments of log-amplitude PSD maximum modeling

To quantify the relevance of the proposed formulations for the log-amplitude PSD maximum, a specific study has been led.

Figure 3 shows the relative error between the frequency of the maximum of the temporal log-amplitude spectrum (computed numerically from Eq. (6)) and the frequency estimated with equations (27a) and (27b), for the polar (Fig. 3a) and equatorial (Fig. 3b) configurations. The relative error is represented versus the main anisotropy ratio $A_z$ and the PSD slope $p_m$, reported in the abscissa and ordinate axes respectively. Let's remark that, in theory, we deal with integral quantities, such as the asymptotic log-amplitude variance, which is valid for $2 < p_m < 6$ (cf. Eq. (24) in Galiègue et al., 2016), or the derivation of the low-frequency asymptote of the log-amplitude PSD, valid for $1 < p_m < 5$ (eq. (13)). These integrals are defined and finite for a restricted range of the spectral index $p_m$. In those cases, the real range of $p_m$ is more restricted than that found in Figure 3. For the polar configuration (left graph), a region of high error is visible for low $A_z$ values (below 5) and $p_m$ values between 3 and 4. In this area, the error reaches nearly 100%, indicating a significant divergence between the frequency of the maximum estimated numerically from (6) and the estimated one with (27a) and (27b). From this area, when $A_z$ or $p_m$ increases, the error gradually decreases. For higher values of $p_m$ and $A_z$, the error becomes very low, indicating a better match between the frequencies. For the equatorial configuration (right graph), the relative errors remain generally lower compared to the polar configuration. The error does not exceed 60%, even in the areas of higher error. The error remains relatively stable as $A_z$ increases, while it gradually decreases with increasing

$p_m$, following a similar trend to the polar configuration but with fewer extreme variations. The region where $p_m$ is close to 5 ($4.5 \le p_m \le 5.5$) gives a high error of 60%, when $A_z \ge 11$ for polar configuration and $A_z \ge 3$ for equatorial configuration. In both configurations, the error is lower for higher values of $A_z$ and $p_m$, suggesting that the frequency estimation is more accurate for these values. On the other hand, the polar configuration appears more sensitive to variations in $A_z$ and $p_m$, showing much higher errors than the equatorial one for low values of $A_z$ and $p_m$. Even if for some combinations of the parameters describing the ionospheric irregularities, the results lead to a non-negligible error, these graphs suggest that the equation applied to estimate the frequency maximum is more reliable than the Fresnel frequency commonly used and could serve in a PSD inversion objective.

A last error study is proposed by considering the relative error on the PSD maximum, versus the outer scale of turbulence $L_0$ and slope $p_m$, reported in abscissa and ordinate axes, respectively. As shown in Figure 4, the estimated frequency using relations (27a) and (27b) provides very good results for low PSD slope ($p_m \le 3.5$) and equatorial configuration, or for high PSD slope ($p_m \ge 5.5$) and polar configuration. However, for slope values $p_m$ close to 4 (for polar configuration) and close to 4.7 (for equatorial configuration), the relative error reaches up to 90% for large outer scale values.

Recalling that these results were obtained assuming $A_y = 1$ the tests were conducted for cases where $1 < A_y \ll A_z$ showed similar performance.

## 3 Comparison of the different approaches

In Sections 2.1 and 2.2, different approaches to model the log-amplitude and phase PSD for a wave illuminating the ionosphere are introduced: the incident plane wave (PW), the incident spherical wave (SW) and a heuristic approach consisting of correcting the Fresnel length in the filtering functions: Corrected Plane Wave (CPW). These approaches are compared
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Relative error in percentage between the frequency of the maximum of the log-amplitude spectrum estimated numerically from (6) and the estimated one from equations (27a) and (27b) in color scale for both polar (a) and equatorial (b) configurations according to the outer scale of turbulence $L_0$ and the slope of the irregularities PSD $p_m$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Temporal power spectral densities for different approaches.</figcaption></figure>

hereafter, and the sensitivity of the PSD to ionospheric parameters is also studied.

The comparison exercise allows studying the received signal PSD for an incident wave illuminating irregularities that are a plane wave or a spherical wave. Figures 5a and 5b show the log-amplitude and phase spectra computed from equations (6) and (7) for the three analytic approaches, plane wave (PW) (in blue), corrected plane wave (CPW) (in red and yellow crosses), and spherical wave (SW) (in red and yellow), in polar configuration defined by Table 1. The red and the yellow curves correspond to a LEO satellite at an altitude of 600 km and a MEO satellite at an altitude of 20000 km respectively. The Corrected Plane Wave (CPW) approach gives a good approximation of the spherical one (SW) for both LEO and MEO satellites. The Fresnel frequencies $f_F^{\text{PW}}$, $f_F^{\text{MEO}}$, and $f_F^{\text{LEO}}$ (Eq. (18)) are indicated with vertical black dashed lines on Figure 5a and do not match with the maximum of the log-amplitude

spectrum. The frequencies of the maximum of the log-amplitude spectrum $f_{\text{max}}^{\text{PW}}$, $f_{\text{max}}^{\text{MEO}}$ and $f_{\text{max}}^{\text{LEO}}$ obtained with equations (27a) and (27b) are plotted in colored dashdotted vertical lines and give a more accurate approximation of the frequencies of the log-amplitude maximum PSD than the Fresnel frequencies. However, it can be observed that these values do not correspond perfectly to the maximum.

Figure 5 allows to study the PW approximation by comparison to SW approach. The spectra for incident plane (PW) and incident spherical wave (SW) approaches present differences for a satellite in Low-Earth Orbit. The shape of log-amplitude and phase PSD is similar. However, a frequency and level shift can clearly be observed on both spectra. The difference is clearly observed in the example plotted in Figure 5a. The low frequency asymptote, the maximum value, and the frequency of the maximum are different between PW and LEO SW. Note that the spectra for the incident plane (PW) and the incident
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Log-amplitude PSD of SW approach, for two configurations, polar and equatorial, and different values of anisotropy $A_z$ and slope $p_m$.</figcaption></figure>

spherical wave (SW) approach converge as the satellite's altitude increases, i.e. when the SW and PW Fresnel distances converge. From this finding, it is preferable to use the more general formalism of the incident spherical wave, in particular for satellites in Low-Earth Orbit. This result was expected and has also been demonstrated in the first part of this study (Morel et al., 2025).

To go further in the analysis, we tested other configurations whose results are plotted in Figure 6. Figures 6a and 6b show the results obtained for a polar and equatorial configurations respectively for different values of anisotropy $A_z$. Polar and equatorial simulation parameters are given in Table 1. In the graphs, the maximum of each log-amplitude PSD is indicated by a point in color. Figures 6c and 6d show the results considering the same configurations, but for different values of slope $p_m$. For the four graphs, the satellite velocity and so the ionospheric piercing point velocity is not taken into account. The impact of the satellite velocity would be to shift the PSD toward higher frequencies, without changing the proposed analysis results. First, we notice that for both

configurations (Figs. 6a and 6b), expression (27b), which is
$$f_{\max} = \frac{V_r}{2\pi} \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)}},$$

corresponds better to the frequency of the maximum than the Fresnel frequency, especially for high anisotropy values. Then, in Figures 6c and 6d, we notice that the highest values of the $p_m$ slope ($p_m = 5$ or 7), the more the frequency of the maximum shifts towards the low frequencies. In this case, the frequency of the maximum could be approximated by expression (27a) which is

$$f_{\max} = \frac{V_r}{2\pi} \sqrt{\frac{4K_0^2}{p_m - 4}}.$$

The opposite reasoning is also valid, when $p_m$ decreases, the frequency of the maximum moves towards high frequencies, and the Fresnel frequency will then correspond better to the frequency of the maximum.

The high and low frequency asymptotes are also plotted in dashed lines in Figure 6. We notice an excellent match between the theoretical spectra and the high-frequency asymptotes for frequencies beyond 5 Hz. The analytic curves match the
low-frequency asymptotes for frequencies lower than 0.1 Hz. These asymptotic elements are particularly relevant for extracting certain parameters characterizing the turbulent ionospheric layer. Indeed, we can observe the dependence of the low frequency asymptote to the anisotropy. Likewise, the peak frequency of the log-amplitude spectrum provides a simple dependence on the speed and the high frequency behavior gives direct access to the inertial slope of the spectrum of inhomogeneities. All these key parameters or formulas could be used in a numerical inversion process of the log-amplitude PSD to retrieve the parameters describing irregularities PSD.

## 4 Conclusion

This paper presents an analysis of the log-amplitude and phase PSD and their potential for estimating parameters that characterize the crossed turbulent ionospheric layer. The study is carried out using an analytic formalism under Rytov's hypothesis, assuming a small turbulent layer thickness compared to the link range.

First, the developed integral formulas enable the investigation of errors potentially induced by using the incident plane wave model (PW) instead of the incident spherical wave model (SW) to predict 3D scintillation effects. By defining realistic propagation scenarios, comparisons of PW and SW approaches to predict scintillation effects are proposed. Initially, it is demonstrated that both approaches (PW and SW) yield similar PSD patterns, albeit with some differences in terms of level, low-frequency asymptote, and the maximum log-amplitude PSD position, particularly when the satellite altitude is low. Furthermore, the approximation consists of replacing the Fresnel distance in the PW formalism of PSD with the SW Fresnel distance is studied. Comparisons with the more rigorous spherical incident wave formalism, indicates that this approximation is highly accurate. Ultimately, this study illustrates the variation trends of the PSD as the satellite altitude changes, with typical including LEO, MEO, and GNSS's typical satellite altitudes.

Additionally, the Rytov formalism allows deriving some asymptotic formulations of the log-amplitude and phase PSD at low and high frequencies. Rytov formalism also allows to study the maximum of the log-amplitude PSD, and several analytic formulations have been proposed. The latter yield satisfactory results in specific configurations with high anisotropy and low value of the slope of the irregularities PSD. The relative error compared to the numerically estimated frequency of the maximum of log-amplitude PSD is calculated to be less than 5% in favorable cases.

This expression of the frequency of the log-amplitude PSD maximum, along with the high- and low-frequency PSD asymptotes, constitute a set of analytic equations useful for addressing the inverse problem. This involves evaluating key parameters describing irregularities PSD, such as the anisotropy ratio along the Earth's magnetic field, the slope of the irregularities PSD and the drift speed of the ionosphere.

These advancements lay the groundwork for an inversion method adapted to recover key parameters of ionospheric medium irregularities. They can be applied to GNSS signals affected by ionospheric scintillation, at high or low latitudes. Another potential application is Synthetic Aperture Radar measurements

performed in the L-band (or at lower frequency), which provide images distorted by ionosphere irregularities.

## Acknowledgments

The authors would like to thank the Direction Générale de l'Arme-ment (DGA), the Agence de l'Innovation de Défense (AID), for their financial support. The editor thanks two anonymous reviewers for their assistance in evaluating this paper.

## Funding

Part of this work was funded by the Agence de l'Innovation de Défense (AID).

## Conflicts of interest

The authors declare no Conflict of Interest.

## Data availability statement

No observational or model data were used for or generated by this work.

## References

* Aarons, J. 1982. Global morphology of ionospheric scintillations. Proc IEEE 70 (4): 360-378. https://doi.org/10.1109/PROC.1982.12314.

* Carrano, CS, Groves KM, Rino CL, Doherty PH. 2016. A technique for inferring zonal irregularity drift from single-station GNSS measurements of intensity and phase scintillations. Radio Sci 51 (8): 1263-1277. https://doi.org/10.1002/2015RS005864.

* Conroy, JP, Deshpande K, Kunduri B, Varney RH, Scales W, Zaghloul A. 2021. Ionospheric scintillation data inversion to characterize the structures associated with a series of polar cap patches. Radio Sci 56 (8): 1-18. https://doi.org/10.1029/2020RS007235.

* Fabbro, V, Feral L. 2012. Comparison of 2D and 3D Electromagnetic Approaches to Predict Tropospheric Turbulence Effects in Clear Sky Conditions. IEEE Trans Ante Propag 60 (9): 4398-4407. https://doi.org/10.1109/TAP.2012.2207070.

* Forte, B, Radicella SM. 2002. Problems in data treatment for ionospheric scintillation measurements. Radio Sci 37 (6): 1-5. https://doi.org/10.1029/2001RS002508.

* Galiègue, H. 2015. Modélisation des effets des scintillations ionosphériques sur la propagation des ondes électromagnétiques en bande L aux latitudes polaires, *Theses*, Université Toulouse 3 Paul Sabatier. https://hal.science/tel-01203660.

* Galiègue, H, Féral L, Fabbro V. 2016. Validity of 2-D electromagnetic approaches to estimate log-amplitude and phase variances due to 3-D ionospheric irregularities. J Geophys Res Space Phys 122: 1410-1427. https://doi.org/10.1002/2016JA023233.

* Ishimaru, A. 1972. Temporal frequency spectra of multifrequency waves in turbulent atmosphere. IEEE Trans Ante Propag 20 (1): 10-19. https://doi.org/10.1109/TAP.1972.1140133.

* Livingston, RC, Rino CL, Owen J, Tsunoda RT. 1982. The anisotropy of high-latitude nighttime F region irregularities. J Geophys Res Space Phys 87 (A12): 10519-10526. https://doi.org/10.1029/JA087iA12p10519.

* McCaffrey, AM, Jayachandran PT. 2019. Determination of the refractive contribution to GPS phase "Scintillation". J Geophys Res Space Phys, 124 (2): 1454-1469. https://doi.org/10.1029/2018JA025759.
Morel, G, Fabbro V, Boisot O, Féral L. 2025. Ionospheric scintillation modeling Part 1: Validity of the plane wave approximation for an electromagnetic wave illuminating the ionosphere. *J. Space Weather Space Clim* **15**: 53. https://doi.org/10.1051/swsc/2025051.

Morel, G, Fabbro V, Boisot O, Ferro-Famil L. 2024. SAR images simulation with ionospheric phase screen. In: *2024 International Radar Conference (RADAR)*, pp. 1–5. https://doi.org/10.1109/RADAR58436.2024.10993564.

Rino, CL. 1979. A power law phase screen model for ionospheric scintillation: 1. Weak scatter. *Radio Sci* **14** (6): 1135–1145. https://doi.org/10.1029/RS014i006p01135.

Rino, CL. 2011. *The Theory of Scintillation with Applications in Remote Sensing*. Wiley-IEEE Press, Hoboken, NJ, USA. ISBN 9780470644775. https://doi.org/10.1002/9781118010211.

Rino, CL, Livingston RC, Tsunoda RT, Robinson RM, Vickrey JF, Senior C, Cousins MD, Owen J, Klobuchar JA. 1983. Recent studies of the structure and morphology of auroral zone F region irregularities. *Radio Sci* **18** (6): 1167–1180. https://doi.org/10.1029/RS018i006p01167.

Rytov, SM, Kravtsov YA, Tatarskii VI. 1989. *Principles of statistical radiophysics. 4. Wave propagation through random media*. Springer Berlin, Heidelberg.

Song, K, Meziane K, Kashcheyev A, Jayachandran PT. 2022. Multifrequency observation of high latitude scintillation: a comparison with the phase screen model. *IEEE Trans Geosci Rem Sens* **60**: 1–9. https://doi.org/10.1109/TGRS.2021.3113778.

Sun, AK, Pi X, Rino C, Lee J. 2023. Statistical uncertainty in the frequency dependence of the intensity scintillation index (S4). *Radio Sci* **58** (7): 1–13. https://doi.org/10.1029/2023RS007659.

Tatarskii, VI. 1971. *The effects of the turbulent atmosphere on wave propagation*. Israel Program for Scientific Translations (IPST), Jerusalem/Springfield, VA, pp. 472. ISBN: 978-0-7065-0680-8.

Taylor, GI. 1938. The Spectrum of turbulence. *Proc R Soc Lond A Math Phys Sci* **164** (919): 476–490. https://doi.org/10.1098/rspa.1938.0032.

Vaggu, PR, Deshpande KB, Datta-Barua S, Bust GS, Hampton DL, Rubio AL, Conroy JP. 2023. Morphological and spectral features of ionospheric structures at E- and F-region altitudes over Poker flat analyzed using modeling and observations. *Sensors* **23** (5): 2477. https://doi.org/10.3390/s23052477.

Vasylyev, D, Béniguel Y, Volker W, Kriegel M, Berdermann J. 2022. Modeling of ionospheric scintillation. *J. Space Weather Space Clim* **12**: 22. https://doi.org/10.1051/swsc/2022016.

Wernik, A, Secan J, Fremouw E. 2003. Ionospheric irregularities and scintillation. *Adv Space Res* **31** (4): 971–981. https://doi.org/10.1016/S0273-1177(02)00795-0.

Wheelon, AD. 2003. *Electromagnetic scintillation*, vol. 2. Cambridge University Press. ISBN 9780511534812. https://doi.org/10.1017/CBO9780511534812.

Yeh, KC, Liu C-H. 1982. Radio wave scintillations in the ionosphere. *Proc IEEE* **70** (4): 324–360. https://doi.org/10.1109/PROC.1982.12313.

**Cite this article as:** Morel G, Fabbro V, Boisot O & Féral L. 2025. Ionospheric scintillation modeling II. Theoretical study of power spectral densities of ionospheric scintillation. *J. Space Weather Space Clim.* **15**, 54. https://doi.org/10.1051/swsc/2025050.

# Appendix A: Asmptotic formulations for log-amplitude and phase spectra

## A.1 Log-amplitude and phase spectra high frequency asymptote

For $\omega$ tends to infinity, $\lim_{\omega \to \infty} \mathcal{F}_\chi$, $\mathcal{F}_\varphi = 1$, equations (6) and (7) becomes:

$$ W_{\chi,\varphi}(\omega) = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p} A_y A_z C_s}{V_u} \int_{-\infty}^{+\infty} \left( A \left( \frac{\omega - k_v V_v}{V_u} \right)^2 + B k_v^2 + 2 C k_v \left( \frac{\omega - k_v V_v}{V_u} \right) + \frac{K_0^2}{a_x^2} \right)^{-p_m/2} dk_v \\ = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p} A_y A_z C_s}{V_u} \int_{-\infty}^{+\infty} (P k_v^2 + Q k_v + R)^{-p_m/2} dk_v, $$

with:

$$ P = A \frac{V_u^2}{V_v^2} + B - 2C \frac{V_v}{V_u}, $$

$$ Q = 2A\omega \frac{V_v}{V_u^2} - 2C \frac{\omega}{V_u}, $$

$$ R = A \frac{\omega^2}{V_u^2} + \frac{K_0^2}{a_x^2}. $$

The integral can be calculated if $P k_v^2 + Q k_v + R > 0$ which implies $4PR - Q^2 > 0$, and if $p_m > 1$. By completing the square, by setting $T^2 = \frac{R}{P} - \frac{Q^2}{4P^2}$, and using the change of variable $x = k_v + \frac{Q}{2P}$, the integral becomes:

$$ \int_{-\infty}^{+\infty} (P k_v^2 + Q k_v + R)^{-p_m/2} dk_v = 2P^{-p_m/2} \int_{0}^{+\infty} (x^2 + T^2)^{-p_m/2} dx, \quad (A.1) $$

which is equal to $P^{-p_m/2} T^{1-p_m} \sqrt{\pi} \frac{\Gamma(\frac{p_m-1}{2})}{\Gamma(p_m/2)}$, after some calculations, the expression of equation (11) appears.

## A.2 Log-amplitude spectrum low frequency asymptote

For $\omega$ tends to zero, $K_0$ is neglected, and simplified weighting functions are used, i.e., only the cosine term is considered (we consider the thin layer assumption i.e., $R_{\text{iono}} / 2L_v \ll 1$). For an incident plane wave, the simplified functions can be expressed as: $\mathcal{F}_\chi(\frac{\omega-k_v V_v}{V_u}, k_v) = 1 - \cos(\frac{L_v}{k_0} (\frac{\omega-k_v V_v}{V_u^2} + k_v^2))$ which is equal to $1 - \cos(\frac{V_u^2 L_v}{V_u^2 k_0} k_v^2)$ when $\omega$ tends to zero.

$$ W_{\chi}(\omega \to 0) = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p} A_y A_z C_s}{V_u} $$
$$ \int_{-\infty}^{+\infty} \left( A \frac{V_v^2}{V_u^2} k_v^2 + B k_v^2 - 2 C k_v^2 \frac{V_v}{V_u} \right)^{-\frac{p_m}{2}} (1 - \cos(ak_v^2)) dk_v $$

Factoring $kv$ and substituting $u = k_v^2$, the previous equation becomes:

$$ W_{\chi}(\omega \to 0) = $$

$$ \frac{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s}{V_u \left( A \frac{V_v^2}{V_u^2} + B - 2C \frac{V_v}{V_u} \right)^{p_m/2}} \int_0^{+\infty} u^{-(p_m+1)/2} (1 - \cos(au)) du, $$

and $\int_0^{+\infty} x^{-\mu}(1 - \cos(bx))dx = \frac{h^{\mu-1}\pi}{2\Gamma(\mu)\cos(\pi\frac{a_x}{2})}$, for $1 < \mu < 3$ and $b > 0$ which is equivalent to $1 < p_m < 5$ and $\frac{V_v^2 L_v}{V_u^2 k_0} > 0$. And so, after some calculations, the expression of equation (13) appears. To get the expression for the incident spherical wave model, just substitute the Fresnel radius for an incident plane wave into the Fresnel radius for incident spherical wave.

### A.3 Phase spectrum low frequency asymptote

The low frequency asymptote for the phase spectrum can now be easily expressed with the log-amplitude asymptote.

$$ W_{\phi}(\omega \to 0) = \frac{\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s}{V_u} $$

$$ \begin{aligned} & \int_{-\infty}^{+\infty} \left( A \frac{V_v^2}{V_u^2} k_v^2 + B k_v^2 - 2 C k_v^2 \frac{V_v}{V_u} + \frac{K_0^2}{a_x^2} \right)^{-p_m/2} (2 - F_{\chi}(...)) dk_v \\ & = \frac{2\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s}{V_u} \int_{-\infty}^{+\infty} \left( A \frac{V_v^2}{V_u^2} k_v^2 + B k_v^2 - 2 C k_v^2 \frac{V_v}{V_u} + \frac{K_0^2}{a_x^2} \right)^{-p_m/2} dk_v - W_{\chi}(0) \\ & = \frac{4\pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s K_0^{-p_m}}{V_u a_x^{-p_m}} \int_0^{+\infty} (\xi k_v^2 + 1)^{-p_m/2} dk_v - W_{\chi}(0) \end{aligned} $$

with $\xi = (A \frac{V_v^2}{V_u^2} + B - 2 \frac{V_v}{V_u}) \frac{a_x^2}{K_0^2}$, the integral is equal to $\frac{B(1/2 p_m/2 - 1/2)}{2\sqrt{\xi}}$, with B being the beta function, And so, after some calculations, the expressions of equation (14) appears.

## Appendix B: Maximum of log-amplitude PSD

### B.1 Basic solutions

In the ionosphere, in particular at low latitudes, the irregularities are strongly anisotropic. The idea is to reduce the two-dimensional spectrum given by (19) to a one-dimension spectrum by determining the direction of the projection of the anisotropy on the $(u, v)$ plane. This direction, where the energy is concentrated is given by $\theta = \frac{1}{2} \operatorname{atan}(\frac{2C}{A-B}) + n \frac{\pi}{2}$, with n an integer. Assuming that the anisotropy is predominant on the z axis, $A_x, A_y \ll A_z$, in this case, it can be shown that $\theta \approx -\alpha_z$. Then,

we can set $k_u = k_r \cos(\theta)$ and $k_v = k_r \sin(\theta)$. For the filtering function, the thin layer approximation allows simplifying the expression (valid for PW and SW formalism):

$$ \mathcal{F}_{\chi}(k_r) \approx 1 - \cos \left( \frac{k_r^2 r_F^2}{2\pi} \right). \quad (B.1) $$

The log-amplitude spectrum becomes:

$$ W_{\chi}(k_r) = \pi r_e^2 \lambda^2 R_{\text{iono}} a_x^{3-p_m} A_y A_z C_s \left( G k_r^2 + \frac{K_0^2}{a_x^2} \right)^{-p_m/2} \left( 1 - \cos \left( \frac{k_r^2 r_F^2}{2\pi} \right) \right), \quad (B.2) $$

with $G = A \cos^2(\theta) + B \sin^2(\theta) + C \sin(2\theta)$. Using the previous assumptions on anisotropies i.e $\theta \approx -\alpha_z$, the expression of G reduces to:

$$ G = A_y^2 \cos^2(\psi) \sin^2(\alpha_y - \alpha_z) + (\sin(\gamma) \sin(\psi) + \cos(\gamma)) \times \cos(\psi) \cos(\alpha_y - \alpha_z)^2. \quad (B.3) $$

Then, using the relationship $\alpha_y = \alpha_z + \arccos(\frac{\tan(\psi)}{\tan(\gamma)})$, which implies $0 < |\psi| \le |\gamma|$ (cf Appendix of Galiègue et al., 2016), it can be shown that G reduces to:

$$ G = A_y^2 \cos^2(\psi) - \frac{\sin^2(\psi)}{\sin^2(\gamma)} (A_y^2 \cos^2(\gamma) - 1) \quad (B.4) $$

Some particular cases give interesting values:

$$ G = \begin{cases} 1, & \text{if } \psi = \gamma \text{ or } A_y = 1, \\ A_y^2, & \text{if } \psi = 0. \end{cases} \quad (B.5) $$

It follows from the differentiation of (B.2) that it is proportional to:

$$ \frac{dW_{\chi}}{dk_r}(k_r) \propto k_r \left( Gk_r^2 + \frac{K_0^2}{a_x^2} \right)^{-p_m/2-1} \left[ \frac{r_F^2}{\pi} \sin \left( \frac{k_r^2 r_F^2}{2\pi} \right) \left( Gk_r^2 + \frac{K_0^2}{a_x^2} \right) + Gp_m \left( \cos \left( \frac{k_r^2 r_F^2}{2\pi} \right) - 1 \right) \right], \quad (B.6) $$

which leads to the following equation system:

$$ \left\{ \begin{array}{l} k_r \\ \left( Gk_r^2 + \frac{K_0^2}{a_x^2} \right)^{-p_m/2-1} \\ \frac{r_F^2}{\pi} \sin \left( \frac{k_r^2 r_F^2}{2\pi} \right) \left( Gk_r^2 + \frac{K_0^2}{a_x^2} \right) + Gp_m \left( \cos \left( \frac{k_r^2 r_F^2}{2\pi} \right) - 1 \right) \end{array} \right. = 0, $$

We are therefore interested in the third equation only, since the first is trivial and the second has no solution because $G > 0$. To solve the third equation, we use the trigonometric formula $\sin^2(x) = \frac{1-\cos(2x)}{2}$, set $X = \frac{k_r^2 r_F^2}{2\pi}$ and $\Delta = \frac{K_0^2 r_F^2}{2n a_x^2 G}$. The equation becomes:
<figure><img src="image_8.png" /></figure>

a Relative error in percentage between the numerical solution of equation (B.8) and the approximated solution given by (B.9).

<figure><img src="image_9.png" /></figure>

b Relative error in percentage between the numerical solution of equation (B.8) and the approximated solution given by (B.10).

<figure><img src="image_10.png" /></figure>

c Solution chosen, the red area corresponds to the solution from the expansion of tangent while the blue area corresponds to the asymptotic expansion solution.

**Figure B.1.** Relative error between the numerical solution of equation (B.8) and the two solutions obtained (a) and (b) in color scale in function of $p_m$ and $\Delta$. In (c), solution to be chosen in color scale according to $p_m$ and $\Delta$.

$$ \sin(X)(X + \Delta) - p \sin^2\left(\frac{X}{2}\right) = 0. \qquad (B.7) $$

By successively dividing the equation by $\sin(\frac{X}{2})$ and $\cos(\frac{X}{2})$, it comes:

$$ \tan(Y) = \frac{4}{p}Y + \frac{2\Delta}{p}, \qquad (B.8) $$

with $X = 2Y$. We cannot find an exact solution for this type of equation, but it is possible to either do an asymptotic expansion of the n-th root, the first root corresponding to the global maximum of the spectrum, or do a series expansion of the tangent. The expressions obtained with the first and the second method respectively gives:

$$ Y^{(1)} \approx \frac{\pi}{2} - \frac{p_m}{2(\pi + \Delta)}, \qquad (B.9) $$

$$ Y^{(2)} \approx \frac{2\Delta}{p_m - 4}, \qquad (B.10) $$

and the corresponding wave numbers are:

$$ k_{\max}^{(1)} = \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)}}, \qquad (B.11) $$
<figure><img src="image_11.png" /><figcaption>**Figure B.2.** Relative error in percentage between the real frequency of the maximum of the log-amplitude spectrum and the estimated one from equation (B.14) and (B.15) in color scale for both polar (a) and equatorial (b) configurations according to the anisotropic ratio $A_z$ and the slope of the irregularities PSD $p_m$.</figcaption></figure>

$$k_{\max}^{(2)} = \sqrt{\frac{4K_0^2}{a_x^2 G(p_m - 4)}}, \quad (B.12)$$

In order to choose the best solution, it is necessary to study the relative error induced by each of the two solutions with respect to the numerical solution.

Figures B.1.a and B.1.b show the relative error for each solution $Y^{(1)}$ and $Y^{(2)}$ with respect to the numerical solution of equation (B.8) as a function of $p_m$ and $\Delta$. One can observe that the validity domain of each solution is limited and corresponds to different areas. The asymptotic development of the n-th root provides better results for low values of the slope $p_m$, while for high values, the series expansion of the tangent is more accurate. For slope values $p_m \le 4$, the solution using the series expansion of the tangent has no solution. Then, the validity domains have been compared to research the optimal solution to apply versus $p_m$ and $\Delta$. Figure B.1.c shows the best solution to choose among both to minimize the relative error. To select the best model versus $p_m$ and $\Delta$ values, it is necessary to determine the equation of a straight line $p_m = f(\Delta)$. Numerically, by polynomial approximation of order 1, we find $p_m \approx 4+2.3\Delta$. Thus the final solution is given by:

$$k_{\max} = \begin{cases} \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)}}, & \text{if } p_m > 4 + 2.3\Delta, \\ \sqrt{\frac{4K_0^2}{a_x^2 G(p_m - 4)}}, & \text{otherwise.} \end{cases} \quad (B.13)$$

(a) Relative error in percentage between the numerical solution of equation (B.8) and the approximated solution given by (B.9). (b) Relative error in percentage between the numerical solution of equation (B.8)

and the approximated solution given by (B.10). (c) Solution chosen, the red area corresponds to the solution from the expansion of tangent while the blue area corresponds to the asymptotic expansion solution.

## B.2 Solutions with one Newton-Raphson iteration

By taking equation (B.8) and by applying one iteration of Newton-Raphson algorithm, the solution obtain with the first and the second method, called $k_{\max}^{1NR}$ and $k_{\max}^{2NR}$ respectively are:

$$k_{\max}^{(1NR)} = \sqrt{\frac{4\pi^2}{2r_F^2} - \frac{4\pi p_m}{2r_F^2(\pi + \Delta)} - \frac{4\pi}{r_F^2} \frac{\tan\left(\frac{\pi}{2} - \frac{p_m}{2(\pi + \Delta)}\right) - \frac{4}{p_m}\left(\frac{\pi}{2} - \frac{p_m}{2(\pi + \Delta)}\right) - \frac{2\Delta}{p_m}}{\sec^2\left(\frac{\pi}{2} - \frac{p_m}{2(\pi + \Delta)}\right) - \frac{4}{p_m}}}, \quad (B.14)$$

$$k_{\max}^{(2NR)} = \sqrt{\frac{4K_0^2}{a_x^2 G(p_m - 4)} - \frac{4\pi}{r_F^2} \frac{\tan\left(\frac{2\Delta}{p_m - 4}\right) - \frac{8\Delta}{p_m(p_m - 4)} - \frac{2\Delta}{p_m}}{\sec^2\left(\frac{2\Delta}{p_m - 4}\right) - \frac{4}{p_m}}}. \quad (B.15)$$

The analytic expressions are therefore more complex and less suitable for inversion parameters of the ionospheric layer. Furthermore, Figure B.2 shows the relative error between these solutions and the numerical solution of the maximum of the log-amplitude spectrum in the configuration of Table 1. By comparing Figure B.2 and Figure 3, it is clear that the improvement in accuracy is not significant, and thus this approach will not be employed.