Method of spherical phase screens for the
modeling of propagation of diverging beams
in inhomogeneous media
Michael E. Gorbunov
1,2
, Oksana A. Koval
1,
*
, Victor A. Kulikov
1
,
and Alexey E. Mamontov
1
1
A.M. Obukhov Institute of Atmospheric Physics, Russian Academy of Sciences, Pyzhevsky per., 3,
Moscow 119017, Russia
2
Hydrometeorological Research Center of Russian Federation, Bolshoy Predtechensky per., 11-13,
Moscow 123242, Russia
Abstract. The phase-screen (split-step) method is widely used for the
modeling of wave propagation in inhomogeneous media. Most known is
the method of flat phase screens. An optimized approach based on
cylindrical phase screen was introduced for the 2-D modeling of radio
occultation sounding of the Earth’s atmosphere. In this paper, we propose a
further generalization of this method for the 3-D problem of propagation of
diverging beams. Our generalization is based on spherical phase screens. In
the paraxial approximation, we derive the formula for the vacuum screen-
to-screen propagator. We also derive the expression for the phase thickness
of a thin layer of an isotropic random media. We describe a numerical
implementation of this method and give numerical examples of its
application for the modeling of a diverging laser beam propagating on a 25
km long atmospheric path.
1 Introduction
The method of phase screens has been widely used for the numerical simulation of the
wave propagation of various nature in inhomogeneous media, including the modeling of the
optical (laser) radiation propagation in a turbulent atmosphere [1–5] and the decimeter
waves propagation during radio occultation sounding of the atmosphere [6–8]. This method
is referred to as split-step. This name reflects the fact that the entire inhomogeneous
medium in this method is represented as a sequence of thin layers, and the propagator
describing the propagation of a wave through each layer is approximately written as the
composition of an infinitely thin layer that forms phase distortions of the wave and a
vacuum propagator describing diffraction.
The phase screen method has a fundamental limitation: it does not account for
backscattering. The method of phase screens can also be considered as a finite-dimensional
approximation of the path integral [9]. In problems of modeling of laser radiation in turbulent
media, especially in order to describe the effect of isotropic turbulence, 2-dimensional phase
*
Corresponding author: kov.oksana20@gmail.com
,
(2019)
https://doi.org/10.1051/itmconf /201930
ITM Web of Conferences
30
CriMiCo'2019
1
150
50
2
27
7
© The Authors, published by EDP Sciences. This is an open access article distributed under the terms of
the
Creative
Commons
Attribution
License 4.0 (http://creativecommons.org/licenses/by/4.0/).
screens are used. In modeling of radio occultation experiments, however, significant
optimization of computational costs is achieved by employing the approximation of one-
dimensional phase screens, because in most cases atmospheric inhomogeneities with vertical
scales from hundreds of meters to kilometers are taken into account, while their horizontal
scales significantly exceed the horizontal size of the Fresnel zone.
In the classical version of the method, flat phase screens are used. This leads to
excessive computational costs when describing a diverging wave: the increasing angle
between the screen and the wavefront at the edges of each screen results in oversampling.
In the two-dimensional (2D) modeling of radio occultation experiments, it turned out to be
quite simple to write down a solution for cylindrical 1-dimensional phase screens [7],
which takes into account the shape of the phase front of the incident wave.
In this paper, we generalize this approach and develop the method of spherical phase
screens, using the paraxial approximation. In the Section 2, we derive the basic relations
using the technique of angular spectra [10]. The main result here is the formula for the
propagator in spherical coordinates in the small angle approximation. In Section 3, we show
that the numerical implementation of this method is no more complicated than the case of
flat phase screens, and we give examples of numerical modeling using the model of
isotropic turbulence. In Section 4, we offer our conclusions.
2 Conclusion of basic relations
2.1 Vacuum propagator
As an example, we will use the beam parameters specified for the DELICAT project
(DEmonstration of LIdar based Clear Air Turbulence detection, Demonstration of clear sky
turbulence detection using lidar) [11, 12]: wavelength
354.84 nm, divergence
2
0.3
mrad. The described method can also be used for other beam parameters, within
the framework of the applicability conditions of the approximation used.
We will consider the Cartesian coordinates
,,x y z
and the spherical coordinates
, , ,r
interconnected as follows:
cos cos ,
cos sin ,
sin .
xr
yr
zr
(1)
Let us consider a wave field in space
,,u x y z
with an imposed radiation condition
that implies waves propagating in the direction of the axis
x
. We denote the boundary
condition for the field in the plane
0x
as
0
,u y z
. Hereinafter, we use the Fourier
transform in the following normalization:
, , exp ,
2
, , exp ,
2
k
f f y z ik y ik z dy dz
k
f y z f ik y ik z d d
(2)
where
2/k
is the wave number. Then the wave field in vacuum is written as follows:
0
, , , exp , , ; , ,
2
k
u r u ikS r d d
(3)
,
(2019)
https://doi.org/10.1051/itmconf /201930
ITM Web of Conferences
30
CriMiCo'2019
1
150
50
2
27
7
2
where
, , ; ,Sr
is phase function, having the following form:
22
, , ; , 1 cos cos cos sin sin .S r r
(4)
The maximum absolute values
,
,
, and
are estimated by the value of the half
beam divergence of
0.15×10
−3
. The propagation distance
r
is estimated at 30 km.
We can write the following approximate expression for the phase function:
2 2 2 2
, , ; , 1 .
22
S r r
(5)
Our numerical solution for the field in an inhomogeneous medium will be based on the
split-step method. In the framework of this method, the medium is split in spherical layers,
represented by phase screens centered at the wave source. We will calculate the field on the
next screen using the vacuum propagator and multiplying the resulting field by the factor
that describes the influence of the random medium:
, , exp , , ,
, ; , exp
2
, , exp ,
2
u r r ik r r
kr
d d P r r ikr
kr
u r ikr d d
(6)
where
, , ,rr
is the phase thickness of the phase screen, and
, ; ,P r r
is the
vacuum propagator. In the paraxial approximation, the following expression can be derived:
22
, ; , exp 1 .
2
rr
P r r ik r
r r r r
(7)
This formula follows the standard split-step method of phase screens. Propagation of the
field from screen to screen is performed in the following steps: 1) The spatial spectrum of
the field in the initial phase screen is calculated. 2) The spectrum is multiplied by the
vacuum propagator. 3) The inverse Fourier transform is taken producing the field on the
next phase screen, without taking into account the medium. 4) The field is multiplied by the
phase factor that takes into account the phase perturbation in the medium between the
screens.
The multiplier
exp ik r
can be omitted because it provides a constant phase addition
in each screen. In this case, the vacuum propagator can be written as follows:
22
, ; , exp .
2
r ikr r
P r r
r r r r
(8)
We rewrite the formula (6) in the operator form:
ˆ
, , , , , ,u r r P r r u r
(9)
where
ˆ
,P r r
is the operator describing the screen-to-screen propagation of the field and
possessing the natural group property:
,
(2019)
https://doi.org/10.1051/itmconf /201930
ITM Web of Conferences
30
CriMiCo'2019
1
150
50
2
27
7
3
2 23 1 12 1 13
ˆ ˆ ˆ
, , , .P r r P r r P r r
(10)
It follows from the accuracy of the group property that this propagator is an exact
solution of the approximate form of the parabolic equation, in which
/yr
and
/zr
are
used as transverse coordinates, and the regular decrease in the amplitude of the spherical
wave as
1/ r
is taken into account. An alternative approach to accounting for wave
divergence is possible using lens coordinates [13–15]. Nevertheless, in the referenced
papers, the lens coordinates are used to describe beams in a homogeneous nonlinear
dispersive medium. The possibility of applying this transformation to an arbitrary
heterogeneous medium requires additional research.
The method of spherical phase screens can also be used to study focused and self-
focused beams. Propagator (7) at describes not a diverging, but a converging beam.
2.2 Phase raid in a thin spherical layer
Each layer of a random medium thick
r
between phase screens
r
and
rr
is described
by the function
, , ,rr
, which is the realization of a 2D random field as a function
of
,
for given
r
and
r
. Consider a realization of the 3-D field of the refractive index
,,N x y z
. Then
, , ,rr
is expressed as follows:
, , , cos cos , cos sin , sin .
rr
r
r r N r r r dr
(11)
To calculate this function, we will use the paraxial approximation. Since
2
r
= 0.675
mm, the values of the second order can be neglected. Thus, we arrive at the following
approximate expression:
, , , , , .
rr
r
r r N r r r dr
(12)
We will neglect the regular part of the refractive index, considering it a constant,
producing in each screen a constant phase shift. The field
,,N x y z
will be considered a
3-D statistically homogeneous and isotropic random field with spectral density
NN
κ
, where
,,
x y z
κ
is the three-dimensional vector of spatial
frequencies, and
κ
. For the correlation function, the standard relation holds:
*3
1 2 1 2
exp .
N
N N i d
rr κ r r κ κ
(13)
Expressing explicitly the correlation of the phase path
, , , , , ,r r r r
and inserting
y
r
and
z
r
, we arrive at the expression for the spectral density
of the phase path:
2
2
2
sinc , , .
2
x
N x x
r
r
d
rr
r
μ
(14)
This formula is written in the approximation of a thin layer of a medium with
rr
,
neglecting the beam broadening. For a layer of a medium with a thickness of the order of
,
(2019)
https://doi.org/10.1051/itmconf /201930
ITM Web of Conferences
30
CriMiCo'2019
1
150
50
2
27
7
4
the external scale, the formula can be approximately rewritten like the corresponding
formula for a flat layer [1,2]:
2
2
0, , .
N
r
rr
r
μ
(15)
3 Numerical modeling
3.1 Numerical modeling method
For turbulent fluctuations of the refractive index, we assume the Kolmogorov–Karman
spectrum:
11/6
2 2 2 2 2
0
( ) 0.033 exp / ,
N n m
C
(16)
where
2
n
C
is the structural constant,
00
2/L
,
0
L
is the external scale,
0
5.92 /
m
,
and
0
is the internal turbulence scale.
The numerical algorithm is similar to the case of flat phase screens. Since the vacuum
propagator (8) is written as a multiplier in the Fourier-transformed space, it is numerically
implemented using the fast Fourier transform.
To generate realizations of random uncorrelated phase screens, we use the discrete form
of relation the correlation of the phase path [16]:
2
2
2
22
22
exp ,
44
j l jl
jl
jl
jj ll
ii
NN
N N N N
(17)
where
,
j l j l
and
jl
are discrete Fourier transform of
jl
:
22
exp .
jl j l
jl
jj ll
ii
NN
(18)
Thus,
jl
are uncorrelated random variables with random phases and rms values:
2
2
4
.
jl
jl
NN
(19)
To take into account the fact that the same discretization step according to angular
coordinates for different radiuses corresponds to different spatial scales, we use adaptive
discretization. On each phase screen, the required sampling step is estimated in angular
coordinates, and if the current sampling step exceeds this estimate, the resolution is
doubled. To this end, the Kotelnikov interpolation is used: in the space of Fourier
transforms. The frequency grid is enlarged, keeping the frequency resolution, and
increasing the maximum frequency twice. In the added grid nodes, the Fourier image of the
field is set to 0. After the inverse Fourier transform, we get a field interpolated to the spatial
grid with a half-step discretization.
,
(2019)
https://doi.org/10.1051/itmconf /201930
ITM Web of Conferences
30
CriMiCo'2019
1
150
50
2
27
7
5
