Scintillation properties of a partially coherent vector beam with vortex phase in turbulent atmosphere

We undertake a computational and experimental study of an advanced class of structured beams, partially coherent radially polarized vortex (PCRPV) beams, on propagation through atmospheric turbulence. A computational propagation model is established to simulate this class of beams, and it is used to calculate the average intensity and on-axis scintillation index of PCRPV beams. On comparison with other classes of structured beams, such as partially coherent vortex beams and partially coherent radially polarized beams, it is found that the PCRPV beams, which structure phase, coherence and polarization simultaneously, show marked improvements in atmospheric propagation. The simulation results agree reasonably well with the experimental results. These beams will be useful in free-space optical communications and remote sensing. © 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement


Introduction
In the past several decades, there have been numerous studies on the propagation properties of laser beams in atmospheric turbulence, owing to their important application in remote sensing, laser radar and free-space optical communications (FSOC) [1,2]. On propagation through the atmosphere, the phase of a laser beam will experience random fluctuations due to turbulence, giving rise to intensity fluctuations (scintillation), beam wander, and increased diffractive spreading. The scintillation is one of the primary factors that limits the performance of remote sensing and FSOC, and therefore the ability to reduce scintillation is of utmost importance for the development of such systems.
It has long been known that partially coherent (PC) beams can have reduced turbulence-induced scintillation compared to their fully coherent counterparts [3][4][5][6][7][8][9][10]. The physical mechanism for this increased resistance to atmospheric turbulence can be explained using the coherent mode representation of optical coherence theory [11,12]. In this representation, PC beams can always be expressed as the incoherent superposition of spatially coherent orthogonal modes. A PC beam therefore sends energy via multiple modes through different propagation channels to the receiver. Due to mutual incoherence of these modes, they do not produce an interference pattern at the detector, resulting in a smooth intensity distribution and lower scintillation. Though the earliest studies of PC beams in turbulence focused on beams with a Gaussian degree of coherence, in recent years there have been a number of studies of beams with non-Gaussian shapes of the degree of coherence [13][14][15].
The scintillation can also be affected by other initial beam parameters, such as the amplitude, phase and polarization. In particular, optical fields carrying a vortex phase have been used for high-data-rate transfer in both optical and quantum free-space communications [16,17]. It was demonstrated in [18][19][20][21][22] that vortex beams are also less distorted by atmospheric turbulence. Numerical simulations of vortex beams of topological charge l = +1 in turbulence showed that the scintillation index (SI) is smaller by a factor of 1.5 compared to the case of a conventional Gaussian beam, and decreases further with an increase in topological charge [23].
The use of nonuniform polarization can also significantly improve the performance of beams in turbulence. The effect of polarization on the statistical properties of beams, especially those with a nonuniform state of polarization (known as vector beams), have been extensively investigated [24][25][26][27][28][29][30][31]. In paricular, Gu et al. theoretically studied the scintillation of spatially coherent nonuniformly polarized beams in turbulence [25], and found that such beams can reduce the scintillation through polarization effects alone; these results were confirmed by experiment [26]. Combining the two aforementioned effects, Cheng et al. studied the scintillation of a vector vortex beam, i.e. a radially polarized beam, and found that it has a lower scintillation than a scalar vortex beam under certain conditions [32].
The aforementioned studies considered the atmospheric propagation characteristics of optical beams with nontrivial coherence, phase, and polarization properties independently, or with two of the three quantities structured simultaneously. It is natural to wonder whether a beam which structures all three quantities at the same time will perform even better in the atmosphere. Such partially coherent radially polarized vortex (PCRPV) beams were recently introduced and studied both theoretically and experimentally on propagation through free space and through ABCD optical systems [33]. Those results demonstrated that such beams display peculiar propagation features in contrast to partially coherent radially polarized (PCRP) beams or pure vortex beams. However, because the state of polarization and the vortex structure of a beam depend on its phase, and partial coherence randomizes the phase of the beam, it is not obvious that structuring all three together will produce additional benefits.
Our aim in this paper is to explore the atmospheric propagation characteristics of PCRPV beams, which structure phase, coherence and polarization simultaneously. A computational model is established to take into account both the randomness of the source as well as the atmosphere it propagates through, and the average intensity and SI of the beam in the receiver plane is calculated. Furthermore, we carry out an experimental study of the scintillation characteristics of PCRPV beams passing through turbulence produced by a laboratory hot plate. Our numerical and experimental results indicate that PCRPV beams are less influenced by atmospheric turbulence compared with ordinary PC beams, PCRP beams, and partially coherent vortex (PCV) beams under the same atmospheric conditions, and that beams with all three forms of structuring provide benefits over those with a fewer number of structures.

Theoretical model for PCRPV beams and the computational model for PCRPV beams in the atmosphere
In this section, we first briefly review the theoretical model for the PCRPV beams, and then establish our computational model for simulating such beams in atmospheric turbulence.

Theoretical model for PCRPV beams
Let us consider a quasi-monochromatic, stochastic electromagnetic beam-like field propagating along the z-axis. The second-order statistical properties of the beam at z = 0 (source plane) can be characterized in the space-frequency domain by a 2 × 2 cross-spectral density (CSD) matrix [34], where r 1 = (x 1 , y 1 ) and r 2 = (x 2 , y 2 ) are two independent position vectors in the source plane, and the elements of the matrix are The quantities E 0x and E 0y in Eq. (2) denote two mutually orthogonal components of the random electric field vector along the x and y directions, respectively, and perpendicular to the z-axis. The asterisk and angular brackets respectively represent the complex conjugate and the ensemble average over the source fluctuations. For PCRPV sources with Gaussian Schell-model correlations, the elements can be expressed as where E 0x and E 0y are the deterministic components of the beam field, given by and r j = x 2 j + y 2 j and ϕ j = arctan y j /x j , ( j = 1, 2) are the radial distance and azimuthal angle in a polar coordinate system. Here l is the topological charge of the vortex phase, and w 0 denotes the beam waist size.
The function µ αβ (r 2 − r 1 ) represents the spectral degree of coherence (SDOC) matrix, a measure of the strength of correlations between the two orthogonal components of the field at two points r 1 and r 2 . We consider the Schell-model case, in which the SDOC is only the function of the separation of two points, and assume a specific Gaussian form for each component, with ∆r = r 2 − r 1 , and δ xx , δ yy and δ xy δ yx are the r.m.s widths of the correlation functions of the different components of the field. We further take B xx = B yy = 1; B xy is the complex correlation coefficient between the x and y components of the electric field, and |B xy | = |B yx |. As was discussed in [35], the realizability condition of a PCRP beam is that |B xy | = |B yx | = 1, δ xx = δ yy = δ xy = δ yx ≡ δ 0 , which indicates that the expressions for the four elements of SDOC of the PCRP beams are identical, and one can show that the above realizability conditions are also applicable for the PCRPV beams. By combining Eqs. (3)-(5), the expressions for the four elements of PCRPV beams can be obtained.

Computational approach for numerically simulating the PCRPV sources
Though the cross-spectral density matrix satisfies a pair of wave equations and can in principle be propagated computationally by Fresnel propagation, it is a function of two vector spatial variables, r 1 and r 2 , or four variables total, and directly propagating this matrix is computationally intensive. Instead, we introduce an ensemble of random electromagnetic fields whose average produces the desired CSD, and propagate this ensemble of fields through the system.
There are two methods for numerically creating the ensemble of a partially coherent source, which we refer to as the phase screen (PS) and complex screen (CS) methods [36]. The former employs random phase screens to generate the statistics of the source, while the latter approach uses random complex screens in which both the amplitude and the phase of the source is modified on transmission. The PS method is most convenient for experimental implementation through the use of phase-only spatial light modulators (SLMs). However, only a source with Gaussian-Schell model (GSM) correlations can be generated or rigorously simulated by this approach. The CS method does not suffer this shortcoming, which makes it ideal for numerical simulation, although it is more difficult to implement experimentally.
To synthesize a PCRPV source using the CS method, a single realization of the x or y component of the electric field in the source plane is taken as where T x (r) and T y (r) are the random complex transmittance functions (screens) which determine the characteristics of the spatial coherence of x and y components. E 0α is the incident electric field before the screens. Therefore, the ensemble average over of the first and second-order Equation (7) indicates that the first-order moment of the CS must be zero and that the second-order moment must be equal to the corresponding SDOC of the source.
To generate a single realization of a CS, we first derive the spatial frequency spectra of the SDOC functions, which are given by their Fourier transforms, where f is a spatial frequency vector. On substituting Eq. (5) into Eq. (8) and applying the realizability condition, the elements of the spatial frequency spectra are of the form In the Fourier domain, the CS can be represented as the square of the power spectrum function multiplied by a zero-mean, unit-variance, complex Gaussian random function where R x and R y stand for the complex Gaussian random function. The function T α (f) is now the Fourier transform of T α (r), and represents one realization of the field that produces the SDOC of Eq. (7). In numerical simulations using MATLAB, these random functions can be written as R where N is the grid number and randn is the built-in function returning a pseudorandom, scalar value drawn from a normal distribution with mean 0 and standard deviation 1. Due to the fact that µ xy and µ yx are equal to µ xx and µ yy , the two components of the random function R x and R y in one realization should be the same in the simulation. By performing the inverse Fourier transform of Eq. (10), the CS is finally synthesized. Part I of Fig. 1 illustrates the square modulus of three typical CSs through numerical calculation. It should be noted that other kinds of partially coherent vector sources with correlation functions whose SDOCs are not Gaussian can also be simulated using the CS method, by replacing the spatial frequency spectra shown in Eq. (9) with other functions. It is also to be noted that this CS method is similar to that used to generate phase screens representing the effect of the atmosphere, to be discussed next.

Computational propagation model for the PCRPV beams in atmosphere
One of the most widely used techniques for numerically simulating the propagation of light waves through atmospheric turbulence is the multiple phase screen method [37], which we briefly outline here. In this method, the turbulence is modeled as a collection of thin random phase screens with the desired turbulence statistics. Phase screens are placed along the propagation path at equal intervals ∆z=z/N T , as shown in part II of Fig. 1. ∆z is the distance between two adjacent phase screens, z is the overall propagation distance, and N T is the total number of screens.
The scenario for simulating the propagation of PCRPV beams in atmospheric turbulence is as follows. The incident electric field of the source E (1) x(y) (r 1 ) = E 0x(y) (r 1 )T α (r 1 ) propagates a distance ∆z in free space, and arrives at the first phase screen. Then, the electric field in the plane of the phase screen is modified by a random phase that represents the accumulated turbulence effect over the distance ∆z, E (2) x(y) (r 2 ) =Ē (1) 0x(y) (r 2 ) exp[iθ 1 (r 2 )], whereĒ (1) 0x(y) is the incident electric field just before the phase screen; θ 1 is the accumulated random phase introduced by the turbulence over the propagation distance ∆z. The propagation steps which follow just repeat the first propagation step, until the beam reaches the last phase screen. The incident electric field is expressed as Finally, the light field is focused by a collecting lens and arrives at the detector in the receiver plane. As noted, the method for synthesizing θ n is similar to that for the synthesis of CS described above but uses the power spectrum Φ n (κ) of the turbulence-induced refractive index fluctuations. The relationship between the power spectrum of the turbulence and the spectrum Φ θ (κ) of the phase screen induced by turbulence is given by the following formula, where k = 2π/λ denotes the wavenumber of a light wave with λ being the wavelength; κ ≡ 2π( f x , f y ) is the spatial frequency vector. Following the same procedure for synthesis of the CS for partially coherent beams, Eq. (11) is first multiplied by a random function R n (f n ); we take a Fourier transform of the result. The real or imaginary parts of the result each represent a valid realization of the turbulence screen phase, e.g. θ n (r n+1 ) = Re{F[R n (f n )Φ θ (f n )]}, where Re represents the real part and F is the Fourier transform operation.
In practice, it is important to take into account three characteristic times τ s , τ t and τ d associated with the system. The time τ s is the characteristic time of amplitude and phase fluctuations of the source; τ t represents the characteristic time of phase fluctuations induced by turbulence, and τ d is the integration time of the detector. The relationship between these three characteristic times plays a significant role in the logic of the numerical simulation. In turbulence theory, it is commonly assumed that the characteristic time of source fluctuations is much smaller than the integration time and the time induced by turbulent phase fluctuations, and also that τ t is greater than τ d , i.e., τ s τ d τ t ; this implies that the detector is only sensitive to the fluctuations induced by turbulence, and averages over the fluctuations of the source.
Under these circumstances, the procedure for simulating for the PCRPV beams in turbulence is as follows: (1) The phase screens of turbulence are synthesized and "frozen" in the propagation link. (2) A number K 1 of CS representing the PCRPV source are synthesized, and each of these source realizations are propagated through the simulated turbulence to the detector plane. (3) The K 1 realizations of the intensity are averaged and this average acts as one realization of the PCRPV beam in turbulence. (4) We loop over steps (1)-(3) K 2 times and finally acquire the average statistics of the PCRPV beams in turbulence. In summary, we must perform two averages to get the physical results: one average over the source fluctuations, and one average over the turbulence flutuations. Figure 2 shows the experimental setup for the generation of PCRPV beams and the measurement of such beams after propagating through thermally-induced turbulence; we divide our description into a discussion of source generation and then a discussion of propagation and detection. A Gaussian beam generated from He-Ne laser (λ = 632.8 nm) passes through a beam expander (BE), a reflecting mirror (RM) and a linear polarizer (LP) whose transmission axis is along the y-direction, and is then focused onto a rotating ground glass disk (RGGD) by a lens L1 with focal lens f 1 = 150 mm, producing partially coherent light with Gaussian correlations. In the experiment, the spatial coherence width of the light is controlled through modulation of the beam spot size on the RGGD. Past the RGGD, the partially coherent light is first collimated by a lens L2, and then passes through a Gaussian amplitude filter (GAF), radially polarized converter (RPC) and a spiral phase plate (SPP) with topological charge l = 1 or l = 2, which transforms the partially coherent linearly polarized light beam into a PCRPV beam. Alignment of the SPP is done when the incident field is fully coherent. If we remove the RPC, the generated source becomes a linearly polarized Gaussian Schell-model vortex source, i.e., a PCV beam. If we remove the RPC and SPP both, the generated source becomes a Guassian-Schell model (GSM) beam. To compare with simulation, the "source plane" is taken to lie right after the SPP.

Experimental generation of PCRPV beams and measurement of their statistical properties through turbulence
The generated source is split into two components by a beam splitter (BS). The reflected part passes through a thin lens L4 with focal length f 4 , then arrives at CCD2 which is used to measure the spatial coherence width of the generated source. Both the distance from the source plane to L4 and from L4 to CCD2 is 2 f 4 (i.e., it is a 4 f imaging system). In this case, the spatial coherence width in the image plane of the CCD2 is the same as that in the source plane. The procedure for measuring the spatial coherence width and the degree of coherence can be found in [38]. The transmitted beam passes over an 35 cm × 50 cm electric hot plate and through a collecting lens L3 with focal length f 3 = 400 mm, then arrives at CCD1 (pixel size 4.4 µm × 4.4 µm) located in the focal plane of L3. The CCD1 records the instantaneous intensity distribution and the output signal is sent to a computer to calculate the average intensity distribution and on-axis scintillation index. The electric hot plate produces turbulence through convection, and the strength of the turbulence is controlled by the temperature T of the hot plate. The beam center is about 2.5 cm above the hot plate surface, and the distance from the SPP to L3 is about 1.5 m. In our experiment, the characteristic time τ s of the partially coherent source fluctuations is primarily determined by the rotation speed of the RGGD, and is in the range from 5.0-10.0µs measured by the digital correlator when the spatial coherence width of the source is from 0.1 mm to 0.4 mm. The characteristic time τ t of turbulence-induced phase fluctuations is 50ms when the temperature of hot plate is 160 • C, and the integration time of the CCD1 is set as 5.0ms. Thus, the detector satisfies the "slow" detector condition, i.e., τ s τ d τ t .

Numerical simulation and experimental results
The numerical simulation is taken to match the propagation stage of the generated beam from the source plane to the observation plane (CCD1 plane) shown in Fig. 2. Ten (N T = 10) phase screens are distributed between the source plane and the plane of L3 to simulate the atmospheric turbulence. The physical size of each phase screen is 40 mm × 40 mm and the sample number is 256×256. The power spectrum density of the turbulence used in the simulation is the Kolmogorov spectrum, given by Φ n (κ) = 0.033C 2 n κ −11/3 , where C 2 n is the refractive-index structure parameter. From the plane of L3 to the observation plane, we assume that the turbulence is absent and a two-step propagation algorithm is used to compute the wave field in the observation plane. Thus, a scaling parameter is introduced and the grid size of the observation screen is set as 4.4µm×4.4µm to match the pixel size of the CCD1 used in the experiment. The instantaneous intensity of partially coherent beams is constructed in the observation plane by taking an average over K 1 =400 screens under a frozen realization of turbulence. The number of turbulence realizations is taken to be K 2 =600. The other parameters such as beam width, coherence width and propagation length used in the simulation are same as those in the experiment.
In order to demonstrate the validity of our computational approach, we first determine the intensity distribution of a GSM beam, a PCV beam with l=1 and a PCRPV beam with l=1 in the plane of CCD1 both numerically and experimentally. Figures 3(a)-3(f) are experimental and numerical density plots of the intensity distribution of the three type beams at the detector in free space (C 2 n =0). For comparison, the corresponding cross-lines of three beams at ρ y =0 are also plotted in Figs. 3(g)-3(i). The beam waist size and the transverse coherence width in the source plane measured in the experiment are 2.0 mm and 0.38 mm, respectively. One finds from Figs. 3(a)-3(f) that the average intensity of the three types of beams all display Gaussian-like profiles in the observation plane (far field). As we know, the transverse profile of the PCV beam or the PCRPV beam near the source plane is a doughnut shape due to the phase or polarization singularity in its center, while in the far field the beam shape degenerates to a Gaussian-like profile due to the random fluctuations of the source. As shown in Figs. 3(g)-3(i), the simulation results of three types of the beams agree reasonable well with the experimental results, implying that our computational model is valid in the case of free-space propagation. Figure 4 illustrates the experimental and corresponding simulation results for the average intensity distribution for the three beam classes in the detector plane after passing through turbulence. It can be seen that the average intensities of the type beams are still Gaussian-like profiles in the observation plane, as expected. The beam widths, defined as the FWHM of the profile, are 0.21 mm, 0.23 mm, 0.21 mm for the GSM beam, PCV beam and PCRPV beam with l=1, respectively, compared to 0.17 mm, 0.19 mm, and 0.18 mm in the absence of turbulence (free space). As expected, the beam widths have increased in the presence of turbulence. Comparing the cross-lines of the average intensity for the experimental and simulation results, as shown in Figs. 4(g)-4(i), we also see that there is good agreement. We now turn to the study of the scintillation properties of PCRPV beams in thermally-induced turbulence. The scintillation of light beams is characterized by the scintillation index (SI), defined as the normalized variance of the intensity fluctuations [1] σ 2 (x, y) = I 2 (x, y) / I(x, y) 2 − 1, where I(x, y) is the instantaneous intensity distribution in the observation plane. The angle brackets denote the ensemble average over the realizations of intensity distribution induced by turbulence. In the experiment, the CCD1 first recorded 1000 frames of instantaneous intensity distribution and stored them in the computer memory. Then, all the frames were imported into MATLAB for data processing as follows: each frame represents one realization of the intensity distribution, and each realization can be represented as a matrix I n (x, y), where x and y are the pixel spatial coordinates, and n denotes the realization number, ranging from 1 to 1000. The total intensity distribution I t (x, y) is obtained by just summing over the 1000 realizations, and dividing by 1000. We assume that the on-axis point is the centroid point of the average intensity, which is calculated using the following formulā With the help of MATLAB, the coordinate of the on-axis point (x c ,ȳ c ) can be determined, and the on-axis SI in the receiver plane is acquired via the following formula, where N denotes the number of total frames. Figure 5(a) presents the experimental results of the on-axis SI versus initial coherence width for a GSM beam, a PCV beam (scalar partially coherent vortex beam) with l=1, and a PCRP beam (partially coherent radially polarized beam, no vortex). One can see that the on-axis SI decreases as the coherence width decreases, agreeing with the general observation that partially coherent beams with smaller coherence width are less affected by turbulence. When the coherence width is relatively large (δ 0 > 0.12 mm), the phase or polarization state has a strong impact on the axial SI. For the same coherence width, the PCRP beam has the lowest value of on-axis SI among the three type beams. For example, the value of on-axis SI is 0.010, 0.080 and 0.068 for the GSM beam, the PCV beam and the PCRP beam when δ 0 = 0.38 mm, indicating that the PCRP beam reduces the SI 30.9% relative to the GSM beam and 15.3% relative to the PCV beam. This advantage is lessened as the coherence width decreases.
There is almost no difference in the on-axis SI between the three type beams when the coherence width is smaller than 0.12 mm. This observation can be explained by the fact that the spatial coherence becomes the dominant factor in determining the on-axis SI when the coherence width is small. The experimental results of the on-axis SI of the PCRPV beam with l=1 and l=2 as a function of coherence width are shown in Fig. 5(b). For comparison, the on-axis SI of the PCRP beam (l=0) is also plotted in Fig. 5(b). One finds that for a relatively large coherence width, the PCRPV beam can further reduce the turbulence-induced scintillation, compared to the PCRP beam, which means that the PCRPV beam has an advantage over the PCRP beam for reducing turbulence-induced scintillation. The advantage is further enhanced with the increase of the topological charge.
The corresponding simulation results of the on-axis SI of the GSM beam, PCV beam and PCRPV beam with l = 0, 1, 2 versus the initial coherence width are illustrated in Figs. 5(c) and 5(d). Compared to Figs. 5(a) and 5(b), the numerical results for the variation of the on-axis SI of these beams with the coherence width are consistent with those of the experimental results, which demonstrate that the PCRPV beams indeed are less affected by the turbulence than a PCRP beam. However, in the simulation results, when the coherence width is smaller than 0.22 mm, the value of the SI is almost the same for the the GSM beam and the PCV beam [see in Fig. 5(c)], and almost the same for a PCRPV beam with l = 0, 1, 2 [see in Fig. 5(d)]. Experimentally, the threshold at which the beams had indistinguishable SI was 0.12 mm. This discrepancy between the experimental and theoretical results may be caused by the theoretical turbulence model not agreeing with the actual experimental turbulence. In the experiment, the exact model for the power density spectrum is unknown, and we adopted the widely used Kolmogorov spectrum for the simulations. Furthermore, we modeled the turbulence as a set of discrete phase screens in the simulations, while in actual turbulence it is continuously distributed along the propagation path. However, the simulation results agree reasonably well with the experimental results for larger values of coherence width.

Conclusion
PCRPV beams combine three distinct turbulence-resistance beam features -partial coherence, nonuniform polarization, and vortex structure -into a single beam. Our results indicate that they provide an efficient way to reduce/overcome the detrimental effects of atmospheric turbulence through simultaneous modulation of fluctuations, polarization, and phase. These beams are relatively straightforward to implement, and will be potentially advantageous in laser radar, free-space optical communications and remote sensing.
To study such beams, we introduced a computational propagation model which can account for both the random vector vortex beam as well as the random atmospheric turbulence. A multi-complex screen approach and a multi-phase screen approach were combined to model the partially coherent vector source and atmospheric turbulence, respectively. The average intensity and scintillation of GSM beams, PCV beams and PCRPV beams on propagation through free space and turbulence was studied both through numerical simulation and experiment. It was shown that the simulation results are reasonably consistent with the experimental results, which demonstrates the validity of the computational propagation model. The PCRPV beams showed a consistent advantage over PCRP beams for reducing turbulence-induced scintillation.
Our results indicate that there is much potential in this three-fold structuring of beams for use in atmospheric applications, and that the three types of structuring complement each other to provide the benefits of all three. The PCRPV beams presented here are just one possible class of beams that have this structuring, and our simulation method can be extended to explore the properties of beams with other configurations of coherence, polarization, and phase.