## Abstract

The study of light propagation though the atmosphere is crucial in different areas such as astronomy, free-space communications, remote sensing, etc. Since outdoors experiments are expensive and difficult to reproduce it is important to develop realistic numerical and experimental simulations. It has been demonstrated that spatial light modulators (SLMs) are well-suited for simulating different turbulent conditions in the laboratory. Here, we present a programmable experimental setup based on liquid crystal SLMs for simulation and analysis of the beam propagation through weak turbulent atmosphere. The simulator allows changing the propagation distances and atmospheric conditions without the need of moving optical elements. Its performance is tested for Gaussian and vortex beams.

© 2016 Optical Society of America

## 1. Introduction

Atmospheric turbulence is a major concern in the performance of optical systems such as ground telescopes, remote optical sensing (lidar) and free-space optical communications (FSOC). The disturbance of a beam propagating in open air depends on multiple factors such as the turbulence strength, propagation distance, light wavelength, kind of the beam and its coherence state [1]. The main effects of the atmospheric turbulence are random intensity fluctuations (scintillation), beam wander and beam spreading that significantly degrade the original signal characteristics. The real-life experimental study of beam propagation for different atmospheric conditions and propagation distances is difficult to perform and therefore the study is mostly reduced to theoretical modeling and numerical simulations.

Recently proposed FSOC protocols based on orbital angular momentum (OAM) information encoding using vortex beams [2–4] demand a realistic experimental analysis of their propagation at different distances under different turbulence conditions. It is well-known that the effects caused by weak atmospheric turbulence can be simulated by means of a number of thin random phase screens located at certain positions in the propagation path, see for example [5, 6]. Rotating phase plates [7–9] can be applied in such studies, however, they do not allow in-situ modification of the atmospheric conditions and the propagation distance. In contrast, this can be done applying spatial light modulators (SLMs) [10–12].

Several experimental configurations for simulating atmospheric turbulence using SLMs have been reported in the last decade. For instance, in [13–19] one phase screen is used for simulating satellite optical communication channels and astronomical imaging in presence of turbulence, while in [4, 20] two phase screens are implemented for simulating more complex situations including ground FSOC channels.

In this work we present an optical system based on two phase screens addressed into two SLMs that allow the simulation of the beam propagation under different weak turbulence conditions and propagation distances. In contrast to previously reported systems, the proposed simulator is robust because the SLMs and digital camera are fixed and the variation not only of the atmospheric conditions but also of the FSOC channel length is achieved in a programmable way, thus making the system flexible and practical to mount (the total optical path between the input and output planes of the setup is ∼ 1.5 m). Note that the cost of the system can be reduced by just implementing the two phase screens into a single SLM display, however, a proper beam redirection is required. Here, we consider two synchronized SLMs in order to simplify its design and to avoid complex beam redirection. While, in principle, the proposed setup is similar to one reported in [4], the authors have not exploited its potential considering only one propagation distance and fixed atmospheric strength.

Another goal of this work is the analysis of the distortion of Gaussian and vortex beams by using the proposed atmospheric turbulence simulator in the context of free-space communication systems. In particular, the intensity fluctuations and spreading of the beams are studied. The choice of a Gaussian beam is twofold. First of all, the theory of its propagation through turbulent atmosphere is well established [5] and will be used for comparison with the experimental results. The scintillation index (SI) experimentally obtained is analyzed not only in the beam center, but also along the transverse beam coordinate. The second reason is that one can expect that the correct experimental simulation of Gaussian beam propagation can be extended to other types of Gaussian-like beams such as Laguerre-Gaussian (LG) or their combinations, which are applied in recently proposed FSOC protocols [2, 4, 9, 21]. As an example, the SI of a vortex LG beam is considered.

## 2. Theoretical background

The refractive index variations, caused by temperature gradients, wind, etc. yield to beam distortions during its propagation through the atmosphere. The strength of the turbulence is often characterized by the refractive-index structure parameter
${C}_{n}^{2}$ which takes values in the range 10^{−17} − 10^{−12} m^{−2}^{/}^{3}, where the higher value yields stronger fluctuations. The statistical behavior of the turbulence is described by the power spectral density of the refractive index fluctuation, Φ(*κ*). There are different models as, for example, the Von Kármán one applied in our study, which retains many physical aspects of the atmospheric behavior. Specifically, the Von Kármán power spectral density model is given by [5]

*κ*= 5.92

_{m}*/l*and

_{i}*κ*

_{0}= 2

*π/l*, with

_{o}*l*and

_{i}*l*being the inner and outer scales of the atmospheric turbulence whose typical values are 0.1cm

_{o}*< l*1cm and 1m

_{i}<*< l*1km. Here

_{o}<*r*

_{0}is an atmospheric coherence width, also known as Fried parameter, defined for the case of plane waves as ${r}_{0}={(0.4229{k}^{2}{C}_{n}^{2}L)}^{-3/5}$, where

*L*is the propagation distance,

*k*= 2

*π/λ*, and

*λ*is the light wavelength. The atmospheric turbulence is classified as weak or strong depending on the value of the Rytov variance ${\sigma}_{R}^{2}=1.23{k}^{7/6}{C}_{n}^{2}{L}^{11/6}$. Note that the weak fluctuation regime corresponds to ${\sigma}_{R}^{2}<1$. In this case, a measurable quantity is the normalized variance of intensity (

*I*) often known as scintillation index (SI): where $\u3008\cdot \u3009$ stands for ensemble average. For the case of plane waves the SI coincides with ${\sigma}_{R}^{2}$. The scintillation is a random intensity variation caused by light refraction and diffraction by turbulent eddies. For Gaussian beams the connection between the SI and ${\sigma}_{R}^{2}$ is more complex: The scintillation index varies along the transverse beam coordinate and its lowest value is at the beam center [5]. Another important parameter is the beam radius

*W*which describes the beam spreading caused by the turbulence in addition to diffraction. This parameter is closely related to the beam wander variance, that measures the root mean squared of beam’s centroid deviation. The theoretical expressions of scintillation and beam spreading considered in this work (see appendix) have been previously reported in [5].

## 3. Experimental simulator design

In the case of weak turbulence, the randomization of the beam phase caused by the the atmosphere’s refractive index variations can be simulated applying phase screens [6]. It has been shown [7] that at least two phase screens are required for proper simulation of beam scintillation and spreading occurring in the medium. Following this approach, an optical communication channel [see Fig. 1(a)] of length *L* and wavelength *λ* operating in a weak turbulent medium characterized by the constant structure parameter
${C}_{n}^{2}$ corresponding to the Fried parameter *r*_{0} can be represented as a pair of thin random phase screens placed at certain distances *L*_{1} and *L*_{2} from the input plane [see Fig. 1(b)]. The Fried parameters *r*_{01} and *r*_{02} of the screens and their positions satisfy the following Eqs.

*z*

_{1}, between phase screens as

*z*

_{2}, and between the latter screen and the receiver as

*z*

_{3}, where

*z*

_{1}+

*z*

_{2}+

*z*

_{3}=

*L*the Fried parameters for the screens are expressed as

Because the parameters *r*_{01} and *r*_{02} have to be positive, there are two constraints for the phase screen positions:

For a concrete propagation distance, screen positions and turbulence strength
${C}_{n}^{2}$, the phase screen distributions for turbulence simulation can be calculated using the programs developed in [6]. In the case when it is preferable to have equal Fried parameters for both screens:
$1.0934{L}^{5/6}={({z}_{2}+{z}_{3})}^{5/6}+{z}_{3}^{5/6}$. If it is desired to have equidistant location of the screens, *z*_{1,2,3} = *z* = *L/*3, then *r*_{01} = 1.5775*r*_{0} and *r*_{02} = 1.4600*r*_{0}.

To setup the laboratory system based on a pair of phase screens [see Fig. 1(c)], it is necessary to perform a proper scaling [4] of the free-space intervals *z _{i}*, beam and screen parameters of the discretized atmospheric model [see Fig. 1(b)]. This can be achieved by considering the

*ABCD*matrix formalism. It is known that the beam propagation in a free-space interval

*z*is equivalent to the propagation of a beam with a smaller size through the laboratory system comprising a free-space interval

*z′*embedded into two spherical lenses as it follows:

*f*and

_{in}*f*are the focal length of the lenses. Here, the parameters

_{out}*α*and

_{in}*α*are transverse coordinate scaling factors defined as the ratio ${\alpha}_{in}={{W}^{\prime}}_{in}/{W}_{in}({\alpha}_{out}={{W}^{\prime}}_{out}/{W}_{out})$ between the input (output) beam radii in the laboratory system ${{W}^{\prime}}_{in}({{W}^{\prime}}_{out})$, and in the real life system

_{out}*W*(

_{in}*W*). From the expression Eq. (7) the following constraints are derived: (8)

_{out}This process is repeated for each of the three free-space intervals *z _{i}* that are associated with the laboratory sub-systems: Input-plane–SLM

_{1}, ${{z}^{\prime}}_{1}$, SLM

_{1}–SLM

_{2}, ${{z}^{\prime}}_{2}$, and SLM

_{2}

*–*output-plane, ${{z}^{\prime}}_{3}$, see Fig. 1(c). In the laboratory system the output scaling at each stage is equal to the input scaling of the next one, i.e:

*α*=

_{out,i}*α*

_{in,i}_{+1}=

*α*, and therefore, the focal length of a lens at the plane

_{i}*i*is given by

*f*=

_{i}*f*

_{in,i}_{+1}

*f*/(

_{out,i}*f*

_{in,i}_{+1}+

*f*). The last lens at the output plane is not needed because only the intensity distribution is measured. Note that the spatial characteristics of the two phase screens implemented by the SLM

_{out,i}_{1}and SLM

_{2}: Fried parameters, the inner and outer scales appearing in Eq. (1) are likewise scaled by its respective factor

*α*

_{1,2}.

Finally the parameters describing the proposed laboratory *simulator* for the beam propagation through a turbulent atmosphere satisfy the system of Eqs.:

*z*are the distances between the phase screens and input/output plane, corresponding to the discretized model of optical communication channel in turbulent atmosphere, as it is displayed in Fig. 1(b), and ${\alpha}_{0}={{W}^{\prime}}_{0}/{W}_{0}$ is the ratio between the input beam radii in the simulator ( ${{W}^{\prime}}_{0}$) and in the real life

_{i}*W*

_{0}. The beam curvature

*F*

_{0}is scaled as ${{F}^{\prime}}_{0}={\alpha}_{0}^{2}{F}_{0}$, correspondingly. Note that in addition the constraint Eq. (6) has to be satisfied.

This system of Eqs. can be simplified for some cases. For instance, one can choose the condition
${{z}^{\prime}}_{1}/{z}_{1}={{z}^{\prime}}_{2}/{z}_{2}={{z}^{\prime}}_{3}/{z}_{3}={L}^{\prime}/L=\gamma $ for the distances in the laboratory simulator [Fig. 1(c)] and the discretized model of atmospheric channel [Fig. 1(b)]. Then the constraint Eq. (6) reduces to
${{z}^{\prime}}_{1}<1.06{{z}^{\prime}}_{2}$, while: *α*_{0} = *α*_{2}, *α*_{1} = *α*_{3} = *γ/α*_{0}, and

Because the value *L′* is fixed, the simulated propagation distance *L* is only controlled by *γ*. Note that the sign of the focal distance of the second lens is always opposite to the sign of the first and the third lenses. If, in addition,
$\gamma ={\alpha}_{0}^{2}$ then *α*_{0} = *α _{i}* and

*f*→ ∞ for any index

_{i}*i*, thus, the laboratory system reduces to a scaled copy of the discretized model of atmospheric channel sketched in Fig. 1(b) where $L={L}^{\prime}/{\alpha}_{0}^{2}$. In the latter setup configuration the propagation distance of the simulated FSOC system can only be changed if the input beam size is also changed thus reducing the potential applications of the simulator.

Another possible simplification of Eq. (11) comes from the choice of the phase screen positions ${{z}^{\prime}}_{i}$ in the laboratory. For instance, if ${{z}^{\prime}}_{i}={z}^{\prime}$ then the focal length at the input plane is

In this case changing of not only the atmospheric conditions but also the beam propagation distance is allowed.

The digital lenses are addressed together with the phase screens into the SLMs that allow creating a flexible programmable system for simulation of different turbulence conditions and propagation distances, where all the optical elements (SLMs and digital camera) are fixed. In general, only two SLMs (SLM_{1} and SLM_{2}) are required for the construction of such smart simulator, since the digital lens ℒ_{0} can be implemented by an electronically tunable (optofluidic) lens. Nevertheless, the use of the SLM_{0} also allows generating complex input beams, for example, optical vortices, which are of interest in FSOC [2].

The proposed optical system is suitable for rapid and programmable variation of both the turbulence conditions and beam propagation distances, however, in practice, its performance is limited by the spatial resolution and temporal response of the SLMs (typically tens of milliseconds for conventional liquid crystal SLMs). The values of the scaling parameters *α*_{0} and *α*_{3} are defined by the real-life beam size, the window size and spatial resolution of the SLM_{0} and the digital camera. While, the scaling parameters *α*_{1} and *α*_{2} depend on the required resolution for a correct implementation of the phase screen and the size of the beam at these screens. Indeed, the scaling parameters *α _{i}* cannot differ from each other too much for appropriate phase screen implementation and beam intensity measurements. In particular, in the case corresponding to Eq. (12) a variation of the scaling by factor 2 implies
${\alpha}_{0}^{2}/2\le \gamma \le 2{\alpha}_{0}^{2}$, that yields to the propagation distance interval given by
${L}^{\prime}/2\le {\alpha}_{0}^{2}L\le 2{L}^{\prime}$. For example, for the following parameters

*W*

_{0}= 0.05 m, ${{W}^{\prime}}_{0}={10}^{-3}\mathrm{m}$ and

*L*′ = 1.5 m, the allowed beam propagation distance is limited in the interval

*L*∈ [1875, 7500] m. On the other side, for correct implementation of the digital lens of diameter

*D*, corresponding to a wavefront phase modulation

*πr*

^{2}

*/λf*, the largest spatial frequency

*D/*4

*λf*has to be bigger than 1/2

*p*, where

*p*is the pixel pitch of the SLM. Therefore, the modulus of the focal distance is restricted by

*| f | > Dp/*2

*λ*. For the typical values:

*p*= 10

*µ*m,

*D*= 10 mm and

*λ*= 0.5

*µ*m digital lenses with

*| f | >*10 cm can be implemented.

## 4. Performance of the simulator and experimental results

The experimental simulator for beam propagation through turbulent atmosphere created in our laboratory comprises two synchronized liquid crystal SLMs (SLM_{1} and SLM_{2}, both Holoeye LCR2500 with pixel size of 19 *µ*m) for encoding of the phase screens (PS_{1} and PS_{2}) and the corresponding lenses ℒ_{1} and ℒ_{2} (with focal lengths *f*_{1} and *f*_{2}), see Fig. 1(c). The lens at the input plane, which coincides in our case with the emitter, referred as ℒ_{0} (focal length *f*_{0}), can be implemented by a low cost electrically tunable fluidic lens in the case of simulating turbulent configurations for different values of *L*. In our case, we have used a SLM (Holoeye PLUTO with pixel size of 8 *µ*m) instead of such a varifocal fluidic lens because this also permits a programmable generation of the input signal, for example, Gaussian and vortex beams. Thus, in our case, the SLM_{0} located at the input plane of the system [emitter TX in Fig. 1(c)] encodes the input beam together with the corresponding digital lens ℒ_{0} in a single phase-only hologram, as reported in [22, 23]. The SLM_{0} has been illuminated by a collimated coherent laser beam (Laser Quantum, Ventus, wavelength *λ* =532 nm and operating at 30 mW output power). In the FSOC simulator, the phase screens are displayed in a square area of 512×512 pixels at the center of the SLM displays with a constant rate of 15 frames per second (fps). The distances between the SLMs and the CCD camera (receiver RX) are
${{z}^{\prime}}_{1}=0.50\phantom{\rule{0.2em}{0ex}}\mathrm{m}$,
${{z}^{\prime}}_{2}=0.57\phantom{\rule{0.2em}{0ex}}\mathrm{m}$ and
${{z}^{\prime}}_{3}=0.50\phantom{\rule{0.2em}{0ex}}\mathrm{m}$ as sketched in Fig. 1(c). The CCD camera (Spiricon SP620U, 12-bit gray-level, pixel size of 4.4 *µ*m and 1620×1200 pixels, 15 fps) located in the receiver plane was synchronized with both SLM_{1} and SLM_{2}. Therefore, every image acquisition corresponds to a single turbulence realization.

For every propagation distance and atmospheric condition, 100 pairs of phase screens (PS1 and PS2) were generated following the Fourier series representation with the addition of subharmonics method described in [6] with *l _{i}* = 1cm and

*l*= 20m. To speed up the generation of the phase screens, the computation was performed by using a Graphics Processing Unit (GPU nVidia GeForce GTX550-Ti) which required only 50 s to generate all the 200 phase screens (about four times faster than using CPU). Note that the GPU allows for parallel computing that significantly speeds up the calculation of the fast Fourier transform used in the phase screen generation algorithm. To validate the experimental results, a numerical simulation of the beam propagation through the setup was performed with the same phase screens and device characteristics (pixel size and dynamic range of the SLMs and the CCD camera) used in the laboratory during the studied tests. The beam propagation through the free space intervals ${{z}^{\prime}}_{1}$, ${{z}^{\prime}}_{2}$ and ${{z}^{\prime}}_{3}$ has been simulated by using the spectrum propagation method reported in [24]. At the distances corresponding to the SLM location the complex field amplitude has been multiplied by the proper phase screen. It has been confirmed that a further increase of the realizations (up to 300 phase screens) in the simulation does not significantly change the beam spreading and scintillation index.

_{o}The beam propagation under weakly turbulent atmosphere conditions, in the range
${C}_{n}^{2}\in [1,5]\times {10}^{-16}{\mathrm{m}}^{-2/3}$, for propagation distances *L* =2:5 km and *L* =4 km have been experimentally studied. The choice of the turbulent atmosphere parameter
${C}_{n}^{2}$ is limited by the condition
${\sigma}_{R}^{2}<0.3$ that has to be satisfied for proper representation of the turbulent continuous medium by means of several phase screens for a given propagation distance *L*, as it is reported elsewhere [25]. The considered parameters of the FSOC simulator and turbulence level are listed in the Table 1.

The considered input Gaussian beam is described by the complex field amplitude as

The beam spreading and SI at the receiver plane of the simulator were analyzed. The beam’s amplitude profiles were obtained by integrating the corresponding 2D average intensity distribution of 100 realizations along orthogonal directions followed by the further mean calculation of their square root and normalization. Using Gaussian fitting for these profiles the beam width, *W* was found. Note, that at *r* = *W* the beam intensity is reduced *e*^{2} times with respect to its maximum value in the center. The SI along the transverse coordinate *r* has been calculated using the Eq. (2). Since the precision in SI estimation mostly depends on the number of realizations *N* as
$1/\sqrt{N}$ in order to improve the viability of the results we have followed a procedure similar to one often used for plane wave SI estimation [26]. We have collected 360 radial intensity profiles for each realization and calculated the intensity variance of the whole ensemble along the time coordinate, thus reaching up to *N* = 36000 equivalent realizations. For the presentation of the SI variation the transverse coordinate has been normalized on the beam radius for comparison of the experimental results with the theoretical expectation.

The optical system was first successfully tested in the case of non-turbulent atmosphere simulation, when only the digital lenses ℒ_{1} and ℒ_{2} were addressed into the SLM_{1}and SLM_{2}, respectively. It has been checked previously using additional imaging system that the width
${{W}^{\prime}}_{0}$ and curvature
${{F}^{\prime}}_{0}$ of the input Gaussian beam were in accordance with the programmed ones. In particular, the considered beam width was
${{W}^{\prime}}_{0}=0.87\phantom{\rule{0.2em}{0ex}}\text{mm}$ that in real free-space coordinates corresponds to
${W}_{0}={{W}^{\prime}}_{0}/{\alpha}_{0}=5.9\phantom{\rule{0.2em}{0ex}}\text{cm}$. While, the considered beam curvature was
${{F}^{\prime}}_{0}=0.88\phantom{\rule{0.2em}{0ex}}\mathrm{m}$ (equivalent to
${F}_{0}={{F}^{\prime}}_{0}/{\alpha}_{0}^{2}=3.9\phantom{\rule{0.2em}{0ex}}\text{km}$ in free space) for *L* = 2.5 km, and
${{F}^{\prime}}_{0}=1.22\phantom{\rule{0.2em}{0ex}}\mathrm{m}$ (*F*_{0} = 5.4 km in free space) for *L* = 4 km. Then the intensity distributions at the receiver plane [RX in Fig. 1(c)] have been measured. In Fig. 2(a) and Fig. 2(b) the experimental amplitude beam profile (blue scatter plot) at the input (TX) and output (RX) plane of the setup are displayed, which have been scaled to the real free-space coordinates by *α*_{0} and *α*_{3}, correspondingly. The measured beam’s output amplitude at the receiver is in good agreement with the theoretical expectation (continuous red line plot), see Fig. 2(b). The SI, see in Fig. 2(c), caused by the refresh of the SLM display is more than one order of magnitude lower than what is observed in the presence of turbulence. Therefore, the intrinsic temporal variations of the SLM display (caused by the fluctuations of the liquid crystal cells) do not significantly affect the measurement of the SI.

The next step was to study the simulator performance in the presence of turbulence. In Fig. 3 the Gaussian beam amplitude profiles and SI at the receiver plane are shown for propagation distance *L* = 2.5 km and different turbulence conditions:
${C}_{n}^{2}=1\times {10}^{-16}{\mathrm{m}}^{-2/3}$ in Fig. 3(a),
${C}_{n}^{2}=2\times {10}^{-16}{\mathrm{m}}^{-2/3}$ in Fig. 3(b), and
${C}_{n}^{2}=5\times {10}^{-16}{\mathrm{m}}^{-2/3}$ in Fig. 3(c). The theoretical curves (red color plots) are compared with the results obtained by means of numerical simulations of the proposed simulator performance (black color plots) and with experimental results (blue color plots).

The beam’s amplitude profiles are displayed in Fig. 3 (first and second columns). The measured beam width and the numerically simulated one are about 8–20% and 4–15% lower than the theoretical beam width calculated using Eq. (16) (see appendix), correspondingly. This discrepancy is also present in the results displayed in Fig. 4 corresponding to the simulation of the propagation distance *L* = 4 km. In all these experimental and numerical simulation examples, the difference in the beam width is larger for the larger
${\sigma}_{R}^{2}$. This fact indicates that further research about the phase screen generation method may be required in order to achieve a realistic atmospheric turbulence simulation [27, 28]. Indeed, it is known that the reported methods for phase screen generation do not accurately reproduce the phenomenon of beam wander [28], which is the main reason of the differences between the measured/simulated and theoretically expected beam radius.

We have observed that the experimental, simulated and theoretical [Eq. (14) in appendix] SI curves almost coincide around the axial point (*r/W* = 0) for all the studied cases, see the third column in Fig. 3 for *L* = 2.5 km and in Fig. 4 for *L* = 4 km. Moreover, the characteristic increase of the SI along the radial coordinate *r* is also present in all these examples. There is a reasonably good agreement between the simulated and experimental SI results. However, the discrepancy between the absolute SI values increases with the increment of *r*, which is again more severe for the large
${\sigma}_{R}^{2}$.

Recently proposed FSOC protocols use optical vortices such as the Laguerre-Gaussian (LG) ones, for information encoding, that makes important to study their distortion during propagation in turbulent atmosphere. In contrast to the Gaussian beam, the coherent LG beam –e.g.:
$E={E}_{0}{r}^{l}\mathrm{exp}\left(-{r}^{2}/{W}_{0}^{2}\right)\mathrm{exp}\left(-\mathrm{i}\pi {r}^{2}/\lambda {F}_{0}\right)\mathrm{exp}\left(-\mathrm{i}l\phi \right)$, with *φ* being the azimuthal angle and *l* the topological charge–has a zero intensity point at the beam center due to its phase singularity. In Fig. 5 it is shown the measured 2D intensity distributions of the previously studied Gaussian beam, Fig. 5(a), and a LG vortex, Fig. 5(b), with topological charge *l* = 3, both corresponding to the average of *N* = 100 realizations (the first column) or individual ones (the rest of the columns) of atmospheric turbulence system with the parameters *L* = 2.5 km and
${C}_{n}^{2}=5\times {10}^{-16}{\mathrm{m}}^{-2/3}$. In Fig. 6(a) the measured and simulated amplitude profiles of this LG beam at the receiver (*W*_{0} = 5.9 cm, *F*_{0} = 3 km) are shown. The corresponding SI curves are displayed in Fig. 6(b). The transverse coordinate was normalized on the effective beam radius
${W}_{\text{LG}}=W\sqrt{1+l}$ estimated from the averaged beam amplitude profiles. From the latter results, we concluded that while the form of the experimental SI curve (blue scatter plot) coincides with the simulated one (black dashed line) their absolute values match only in the region of the intensity maximum. We recall that it was also the case for the Gaussian beam. The difference between the simulated and experimental data is bigger for the smaller average intensity values. It could be explained by the background noise that significantly reduces the SI in the intervals with low intensity distribution. Note that the coordinate of the local SI maximum coincides with the inflection point of the averaged intensity profile. This fact as well as the local SI minimum in the vortex singular point (*r* = 0) is worth to verify in the outdoors experiments.

## 5. Discussion and conclusions

In this work we have presented an experimental programmable setup based on two SLMs for simulating the beam propagation through weak turbulent atmosphere for different distances and atmospheric parameters. Its validity is demonstrated with experimental and simulated results for the propagation of Gaussian and vortex beams. To the best of our knowledge it is the first time that the SI is experimentally analyzed not only at the points of maximum intensity (axial point for a Gaussian beam), but along the entire beam profile. We have found that for both kinds of beams the experimental, simulated and theoretical SI curves have the same behavior. Moreover, in the region of maximum intensity, which also coincides with the minimum beam intensity gradient, the experimental, simulated and theoretically predicted results are in a good quantitative agreement for all the studied atmospheric conditions. However, we have observed a discrepancy in the absolute SI values for other regions of the curves especially the ones where the intensity profiles of the experimental, simulated and theoretical beam profiles do not coincide. It seems that the beam intensity gradient plays an important role in the local intensity scintillation. As we have mentioned before, the beam spreading according to the theory is larger than what we obtain in simulations and experiments that also may provide the discrepancy in SI. Probably, the implementation of more sophisticated methods for the phase screen, which use the Zernike polynomials decomposition might help to overcome this problem. We also underline that the theoretical SI curve has a parabolic form [see Eq. (14)] that is certainly not true for large *r* values. Other possible limitations lie in the characteristics (pixel size and the dynamic range) of the SLM used for the phase screen and lens encoding. The observed differences between the simulated and experimental results indicate that the linear SLM response assumed in the simulations is not completely accurate. However, the further progress in the SLM technology will certainly improve the performance of the proposed simulator. Apart from the traditional optical communication schemes the development of recently proposed FSOC protocols based on information encoding into beam shape or its OAM could benefit from the use of this simulator.

## Appendix

In this section the theoretical formulas used to calculate the turbulent parameters of a Gaussian beam are presented. According to [5], the scintillation index using the finite inner and outer scale Von Kármán power spectral density has the form

${\mathrm{\Theta}}_{0}=1-L/{F}_{0}$,
${\mathrm{\Lambda}}_{0}=2L/k{W}_{0}^{2}$,
$\mathrm{\Theta}={\mathrm{\Theta}}_{0}/\left({\mathrm{\Theta}}_{0}^{2}+{\mathrm{\Lambda}}_{0}^{2}\right)$ and
$\mathrm{\Lambda}={\mathrm{\Lambda}}_{0}/\left({\mathrm{\Theta}}_{0}^{2}+{\mathrm{\Lambda}}_{0}^{2}\right)$ are the free space propagation parameters for the input beam radius *W*_{0} and curvature *F*_{0} (for the case *F*_{0} ≠ *L*). The beam spreading due to turbulence is given by

The beam characteristics given by Eq. (14) and (16) are used for comparison with our experimental results.

## Acknowledgments

Spanish Ministerio de Economía y Competitividad is acknowledged for projects TEC2014-57394-P and FIS2013-46475-C3-1-P. Carolina Rickenstorff gratefully acknowledges the CONACyT postdoctoral grant 234821.

## References and links

**1. **G. Gbur and E. Wolf, “Spreading of partially coherent beams in random media,” J. Opt. Soc. Am. A **19**, 1592–1598 (2002). [CrossRef]

**2. **J. Wang, J. Y. Yang, I. M. Fazal, N. Ahmed, Y. Yan, H. Huang, Y. Ren, Y. Yue, S. Dolinar, M. Tur, and A. E. Willner, “Terabit free-space data transmission employing orbital angular momentum multiplexing,” Nat. Photon. **6**, 488–496 (2012). [CrossRef]

**3. **M. Krenn, R. Fickler, M. Fink, J. Handsteiner, M. Malik, T. Scheidl, R. Ursin, and A. Zeilinger, “Communication with spatially modulated light through turbulent air across Vienna,” New J. Phys. **16**, 113028 (2014). [CrossRef]

**4. **B. Rodenburg, M. Mirhosseini, M. Malik, O. S. Magaña-Loaiza, M. Yanakas, L. Maher, N. K. Steinhoff, G. A. Tyler, and R. W. Boyd, “Simulating thick atmospheric turbulence in the lab with application to orbital angular momentum communication,” New J. Phys. **16**, 033020 (2014). [CrossRef]

**5. **L. C. Andrews and R. L. Phillips, *Laser Beam Propagation through Random Media* (SPIE Press, 2005). [CrossRef]

**6. **J. D. Schmidt, *Numerical Simulation of Optical Wave Propagation with Examples in MATLAB* (SPIE Press, 2010) Chap. 9, pp. 149–184.

**7. **S. V. Mantravadi, T. A. Rhoadarmer, and R. S. Glas, “Simple laboratory system for generating well-controlled atmospheric-like turbulence,” Proc. SPIE **5553**, 290–300 (2004). [CrossRef]

**8. **P. Polynkin, A. Peleg, L. Klein, and T. Rhoadarmer, “Optimized multiemitter beams for free-space optical communications through turbulent atmosphere,” Opt. Lett. **32**, 885–887 (2007). [CrossRef] [PubMed]

**9. **Y. Ren, H. Huang, G. Xie, N. Ahmed, Y. Yan, B. I. Erkmen, N. Chandrasekaram, M. P. J. Lavery, N. K. Steinhoff, M. Tur, S. Dolinar, M. Neifeld, M. J. Padgett, R. W. Boyd, J. A. Shapiro, and A. E. Willner, “Atmospheric turbulence effects on the performance of a free space optical link employing orbital angular momentum multiplexing,” Opt. Lett. **38**, 4062–4065 (2013). [CrossRef] [PubMed]

**10. **G. D. Love, “Wave-front correction and production of Zernike modes with a liquid-crystal spatial light modulator,” Appl. Opt. **36**, 1517–1524 (1997). [CrossRef] [PubMed]

**11. **M. A. A. Neil, M. J. Booth, and T. Wilson, “Dynamic wave-front generation for the characterization and testing of optical systems,” Opt. Lett. **23**, 1849–1851 (1998). [CrossRef]

**12. **M. A. A. Neil, F. Massoumian, R. Juškaitis, and T. Wilson, “Method for the generation of arbitrary complex vector wave fronts,” Opt. Lett. **27**, 1929–1931 (2002). [CrossRef]

**13. **C. Paterson, I. Munro, and J. C. Dainty, “A low cost adaptive optics system using a membrane mirror,” Opt. Express **6**, 175–185 (2000). [CrossRef] [PubMed]

**14. **L. Hu, L. Xuan, Z. Cao, Q. Mu, D. Li, and Y. Liu, “A liquid crystal atmospheric turbulence simulator,” Opt. Express **14**, 11911–11918 (2006). [CrossRef] [PubMed]

**15. **L. Jolissaint, “Optical turbulence generators for testing astronomical adaptive optics systems: A review and designer guide,” PASP **118**, 1205–1224 (2006). [CrossRef]

**16. **C. C. Wilcox, J. R. Andrews, S. R. Restaino, T. Martinez, and S. W. Teare, “Atmospheric turbulence generator using a liquid crystal spatial light modulator,” in “2007 IEEE Aerospace Conference,” (2007), pp. 1–8.

**17. **L. Hu, L. Xuan, D. Li, Z. Cao, Q. Mu, Y. Liu, Z. Peng, and X. Lu, “Real-time liquid-crystal atmosphere turbulence simulator with graphic processing unit,” Opt. Express **17**, 7259–7268 (2009). [CrossRef] [PubMed]

**18. **K. Murphy, D. Burke, N. Devaney, and C. Dainty, “Experimental detection of optical vortices with a Shack-Hartmann wavefront sensor,” Opt. Express **18**, 15448–15460 (2010). [CrossRef] [PubMed]

**19. **B. Rodenburg, M. P. J. Lavery, M. Malik, M. N. O’Sullivan, M. Mirhosseini, D. J. Robertson, M. Padgett, and R. W. Boyd, “Influence of atmospheric turbulence on states of light carrying orbital angular momentum,” Opt. Lett. **37**, 3735–3737 (2012). [CrossRef] [PubMed]

**20. **T. L. Kelly, D. F. Buscher, P. Clark, C. Dunlop, G. Love, R. M. Myers, R. Sharples, and A. Zadrozny, “Dual-conjugate wavefront generation for adaptive optics,” Opt. Express **7**, 368–374 (2000). [CrossRef] [PubMed]

**21. **G. Gibson, J. Courtial, M. Padgett, M. Vasnetsov, V. Pas’ko, S. Barnett, and S. Franke-Arnold, “Free-space information transfer using light beams carrying orbital angular momentum,” Opt. Express **12**, 5448–5456 (2004). [CrossRef] [PubMed]

**22. **J. A. Rodrigo, T. Alieva, A. Cámara, Ó. Martínez-Matos, P. Cheben, and M. L. Calvo, “Characterization of holographically generated beams via phase retrieval based on Wigner distribution projections,” Opt. Express **19**, 6064–6077 (2011). [CrossRef] [PubMed]

**23. **V. Arrizón, U. Ruiz, R. Carrada, and L. A. González, “Pixelated phase computer holograms for the accurate encoding of scalar complex fields,” J. Opt. Soc. Am. A **24**, 3500–3507 (2007). [CrossRef]

**24. **D. Mas, J. Garcia, C. Ferreira, L. M. Bernardo, and F. Marinho, “Fast algorithms for free-space diffraction patterns calculation,” Optics Communications **164**, 233–245 (1999). [CrossRef]

**25. **A. Belmonte, “Feasibility study for the simulation of beam propagation: consideration of coherent lidar performance,” Appl. Opt. **39**, 5426–5445 (2000). [CrossRef]

**26. **W. M. Hughes and R. B. Holmes, “Pupil-plane imager for scintillometry over long horizontal paths,” Appl. Opt. **46**, 7099–7109 (2007). [CrossRef] [PubMed]

**27. **E. M. Johansson and D. T. Gavel, “Simulation of stellar speckle imaging,” Proc. SPIE **2200**, 372–383 (1994). [CrossRef]

**28. **J. Recolons and F. Dios, “Accurate calculation of phase screens for the modelling of laser beam propagation through atmospheric turbulence,” Proc. SPIE **5891**, 1–12 (2005).