Optical Field-to-Field Translation under Atmospheric Turbulence: A Conditional GAN Framework with Embedded Turbulence Parameters

1. Introduction

When a laser beam propagates through the atmosphere, it is significantly influenced by small temperature fluctuations, as well as wind shear and turbulence. These factors lead to irregular distributions of the refractive index in both space and time, resulting in the distortion of the laser wavefront. Consequently, this distortion degrades the coherence of the laser, causes dissipation of laser energy, and significantly impacts the transmission efficiency of the laser [1,2,3,4]. Specifically, atmospheric turbulence can cause several detrimental effects on laser communication:

Beam Spread: Atmospheric turbulence causes the laser beam to spread out as it propagates, leading to a larger spot size at the receiver and reduced signal intensity [1,2].

Beam Wander: The random variations in the refractive index cause the laser beam to wander off its intended path, making it difficult for the receiver to maintain alignment and track the beam accurately [1,2].

Angle-of-Arrival Fluctuations: The angle at which the laser beam arrives at the receiver fluctuates due to turbulence, complicating the task of maintaining optimal alignment and reducing the efficiency of signal reception [1,2].

Scintillation: The intensity of the laser beam fluctuates randomly as it propagates through the turbulent atmosphere, leading to signal fading and increased error rates [1,2].

Phase Distortion: Turbulence induces random phase changes in the laser beam, which can degrade the performance of coherent communication systems by reducing the mixing efficiency and increasing error rates [1,2].

Currently, the random phase screen method, based on wave theory, is commonly used to simulate laser beam propagation in a turbulent atmosphere [5,6]. It models the continuous medium along the laser path as discrete phase screens and treats the optical-field propagation between them as vacuum diffraction. Specifically, before reaching each screen, the optical field undergoes vacuum diffraction and is then modulated by phase screens. Iterating this process in a sequential manner can efficiently simulate the laser beam propagation in turbulent media.

The key steps in this numerical simulation encompass the computation of phase screens and the diffraction of the optical field. For the diffraction calculation, it is typically addressed by the split-step method based on the Fresnel diffraction in angular spectrum form, which has become a well-established approach in numerical diffraction computations [7]. Regarding the phase screen calculation, the power-spectrum inversion algorithm is widely utilized to generate phase screens [8,9,10]. Considering its substantial influence on the simulation results, a great deal of research has been to enhance computational efficiency without compromising accuracy, including sparse spectrum technology introduced by Charnotskii [11,12] , random subharmonic sampling by Paulson [13] and sinc series approximation based method introduced by Cubillos[14]. Additionally, the authors have proposed a phase screen generation technique based on sampling optimization [15,16]. All these methods are aimed at improving the speed and precision of phase screen generation by sampling the refractive index spectrum.

Despite these advancements, several challenges persist. Even though remarkable improvements have been made in numerical simulation efficiency, achieving real-time applicability is still extremely difficult. This is mainly attributed to two key factors. Firstly, each transmission calculation involves multiple iterations of diffraction and modulation processes, which are inherently time-consuming operations. Secondly, given the stochasticity of atmospheric turbulence, obtaining the necessary optical field statistics demands statistical analysis of a large number of samples. Typically, this can only be accomplished through ensemble average. These challenges are deeply rooted in the numerical computation approach, indicating that mere algorithmic improvements are insufficient to significantly enhance the simulation speed.

Fortunately, recent breakthroughs in artificial intelligence, specifically generative AI, present promising solutions for these challenges. In the realm of image generation, generative AI has witnessed remarkable progress. It has evolved from variational autoencoders [17] and generative adversarial networks [18,19] (GAN) to the state-of-the-art diffusion models [20], which have substantially enhanced both image quality and content control. In our study, taking the initial optical field distribution as input and atmospheric turbulence conditions as control parameters, generative AI can generate the atmospheric-turbulence-altered intensity distribution on the received plane.

Among the aforementioned generative AI frameworks, GANs stand out for their rapid output capabilities, generating the intensity distribution of laser transmission through atmospheric turbulence more efficiently than diffusion models, which require multiple denoising steps. By avoiding the blurring issues associated with variational autoencoders, GANs are particularly well-suited for the fast field mapping of laser transmission under a turbulent atmosphere. Thus, here we focus on a GAN-based architecture, designing a specialized neural network structure tailored for simulating laser beam propagation in the turbulent atmosphere. Through adversarial training, the authors have achieved fast and stable simulations for laser beam propagation in atmospheric turbulence.

2. Network Design

2.1. The Laser Beam Model

The goal of the proposed model is to simulate the intensity distribution of a laser beam after propagating through turbulent atmospheres. The input to the model is the optical field distribution at the source plane, and the output is the intensity distribution on the receiving plane, while the transmission distance and the refractive index structure constant are provided as additional control parameters in the bottleneck layer.

The mathematical model for laser beam propagation in turbulent atmospheres is based on the split-step method using Fresnel diffraction in angular spectrum form [6,21]. The propagation process can be expressed as:

$$\begin{array}{cc}\hfill U_n({\mathbf{r}}_n,z_n)& =Q[{\displaystyle \frac{m_{n-1}-1}{m_{n-1}{z}_{n-1}}},{\mathbf{r}}_n]\hfill \\ \hfill & \times \prod _{i=1}^{n-1}\left\{T[z_i,z_{i+1}]F^{-1}[{\mathbf{f}}_i,{\displaystyle \frac{{\mathbf{r}}_{i+1}}{m_i}}]Q[{\displaystyle \frac{4{\pi }^2{z}_i}{m_{i}k}},{\mathbf{f}}_i]F[{\mathbf{r}}_i,{\mathbf{f}}_i]{\displaystyle \frac{1}{m_i}}\right\}\hfill \\ \hfill & \times \left\{Q[{\displaystyle \frac{1-m_1}{z_1}},{\mathbf{r}}_0]T[z_0,z_i]U_0\left({\mathbf{r}}_0,0\right)\right\},\hfill \end{array}$$

(1)

where, $U_0({\mathbf{r}}_1,0)$ is the spatial distribution of the optical field at the source plane. m serves as the scaling parameter that controls the spatial sampling interval on the receiving plane, which helps accommodate the beam width increase caused by turbulence and diffraction effects. Q, F, $F^{-1}$, T are a set of operators. Specifically, Q represents the second-order phase factor, F represents the Fourier Transform, $F^{-1}$ represents the Inverse Fourier Transform, and T represents the phase modulation induced by turbulence. They are defined as follows:

$$Q\left[c,\mathbf{r}\right]\left\{U\left(\mathbf{r}\right)\right\}\equiv e^{{\displaystyle \frac{ikc}{2}}{\left|\mathbf{r}\right|}^2}U\left(\mathbf{r}\right),$$

(2)

$$F[\mathbf{r},\mathbf{f}]\left\{U\left(\mathbf{r}\right)\right\}\equiv \underset{-\infty }{\overset{\infty }{\int }}U\left(\mathbf{r}\right)e^{-i2\pi \mathbf{f}·\mathbf{r}}d\mathbf{r},$$

(3)

$$F^{-1}[\mathbf{f},\mathbf{r}]\left\{U\left(\mathbf{f}\right)\right\}\equiv \underset{-\infty }{\overset{\infty }{\int }}U\left(\mathbf{f}\right)e^{i2\pi \mathbf{f}·\mathbf{r}}d\mathbf{f},$$

(4)

$$T[z_i,z_{i+1}]=e^{-i\varphi \left({\mathbf{r}}_{\mathbf{i}+\mathbf{1}}\right)},$$

(5)

where, $\varphi \left(r_{i+1}\right)$ represents the phase screen that can be calculated via phase screen generation method based on spectral inversion:

$$\begin{array}{cc}\hfill \phi (x,y)& ={\sum }_{f_x}{\sum }_{f_y}\mathbf{h}(f_x,f_y)\sqrt[]{{\mathrm{\Phi }}_{\phi }(f_x,f_y)}\hfill \\ \hfill & \times exp\left[j(f_{x}x+f_{y}y)\right]\mathrm{\Delta }{f}_x\mathrm{\Delta }{f}_y.\hfill \end{array}$$

(6)

In Eq.(6), $\mathbf{h}(f_x,f_y)$ represents complex Hermitian Gaussian random matrix with zero mean and unit variance. ${\mathrm{\Phi }}_{\phi }(f_x,f_y)$ denote the power spectral density function that conforms to the atmospheric turbulence spectrum. $f_x$ and $f_y$ are the spatial frequencies in the x direction and y direction, respectively. The power spectral density function that follows the Kolmogorov spectrum can be expressed as:

$${\mathrm{\Phi }}_{\phi }(f_x,f_y)=0.023r_0^{-5/3}{f}^{-11/3},$$

(7)

where $r_0$ represents the atmospheric coherence length, also known as the Fried parameter [1]. It is related to the refractive index structure constant $C_n^2$ and the propagation distance z in a Kolmogorov spectrum:

$$r_0={\left(0.423C_n^2{k}^{2}z\right)}^{-3/5},$$

(8)

2.2. Proposed GAN Architecture

By analyzing the calculation process of the above-mentioned model, it can be deduced that to complete the calculation of laser beam propagation in turbulent atmosphere, it is essential to feed the source-plane optical-field distribution, transmission distance z, and refractive index structure constant $C_n^2$ into the calculation model as inputs, and the final output is the optical-field distribution on the receiving plane. To simulate this process efficiently, we propose a GAN-based architecture. The generator, based on the U-Net [22] network, takes the source-plane optical field distribution as input and outputs the intensity distribution on the receiving plane. The discriminator, based on the PatchGAN network [23], evaluates the quality of the generated samples. The structure of the U-Net generator is shown in Figure 1.

As depicted in Figure 1, the U-Net network essentially adopts an encoder-decoder architecture. The left-hand side of the network progressively abstracts the input optical-field distribution, extracting deep-level features at the middle layer. Subsequently, through the decoding process on the right-hand side, it reconstructs an output optical-field distribution with the same size as the input. Considering that most statistical analyses of atmospheric turbulence predominantly rely on the intensity distribution, the optical-field distribution is replaced by the intensity distribution.

During the GAN training process, a discriminator is necessary to assess the quality of the generated samples. In this study, PatchGAN is used as a discriminator, as it has been proven effective in similar tasks according to previous research [19,24]. The architecture of the PatchGAN network is illustrated in Figure 2.

3. Training Data and Network Training

The training sample was prepared via the numerical simulation method mentioned above. Based on the network architecture, the inputs for the U-Net network are defined as the initial optical field distribution (comprising the amplitude distribution and phase distribution) and the control condition parameters (including the refractive index structure constant and transmission distance).

The initial optical field distribution follows a fundamental mode Gaussian distribution. The initial phase distribution varies depending on the collimation and focusing state of the beam. To match the dimensionality of the middle layer of the U-Net network, the characteristic parameters of the intensity of the turbulence are expanded in dimension to form a characteristic layer (with a dimension of 16 × 16 × 2). Here, the third dimension represents different control condition parameters. Finally, this characteristic layer is concatenated with the feature layer along the third dimension.

As for the PatchGAN network, its inputs consist of the initial optical field distribution and the intensity distribution of the receiving plane. The intensity distribution output of the U-Net network is labeled as ’fake’, while the intensity distribution obtained through numerical simulation is labeled ’true’.

3.1. Building up the Training Dataset

The parameters of the training sample set are presented as follows:

• The input optical field is a fundamental-mode Gaussian beam with a beam waist radius of 0.1 m and a wavelength of 1.08 $\mu$m. This Gaussian beam can be in a collimated or focused state. In the focused state, the focal length is set equal to the corresponding transmission distance;
• Regarding the turbulence-related parameters, the refractive index structure constant is linearly sampled 100 times from $1\times {10}^{-16}$ to $1\times {10}^{-13}$${\mathrm{m}}^{-2/3}$. The outer scale parameter is 100 m, and the inner scale parameter is 0.01 m. The number of phase screens on the transmission path is set to 20;
• For the receiving-plane and transmission parameters, the physical dimension of the receiving-plane is $1\mathrm{m}\times 1\mathrm{m}$. The transmission distance ranges from 0 to 10000 m, with a linear separation of 200 m intervals.

The above-mentioned parameters are substituted into the numerical model based on the random phase screen method to calculate the intensity distribution of the receiving plane. Considering the stochasticity of the atmospheric turbulence, the calculations are repeated 400 times for each set of parameters, and the average result is used as the label data for the network.

Subsequently, a training data set consisting of 10,000 training samples is created. Each sample includes the source-plane optical field with a dimension of 512 × 512 × 2, the control-condition parameters of 16 × 16 × 2, and the label representing the intensity distribution with a dimension of 512 × 512.

3.2. Loss Function

To align the intensity distribution generated by the network with the results of numerical simulation, a loss term representing the statistical characteristics of the atmospheric turbulence effect is incorporated into the original loss function of the generator network. Given that the loss function should be as simple as possible to ensure the gradient stability during backpropagation, a first-order statistical moment, specifically the beam wander variance [1], is introduced as an additional loss term. This additional loss term can be expressed by the following formula:

$$〈r_c〉=\sqrt[]{〈r_{cx}^2〉+〈r_{cy}^2〉}=\sqrt[]{{\left({\displaystyle \frac{\int \int x〈I(x,y,L)〉dxdy}{\int \int 〈I(x,y,L)〉dxdy}}\right)}^2+{\left({\displaystyle \frac{\int \int y〈I(x,y,L)〉dxdy}{\int \int 〈I(x,y,L)〉dxdy}}\right)}^2}.$$

(9)

The total power of the beam in the receiving plane is the sum of the power weighted in the x direction and y direction respectively. When the two-dimensional intensity distribution is obtained, the beam wander variance can be obtained by numerical solution, and then the $L_1$ loss of the beam wander variance can be introduced,

$$\begin{array}{cc}\hfill L_{BW}& =\sqrt[]{{\left({\displaystyle \frac{\sum _{i=-256}^{255}\sum _{j=-256}^{255}i·B_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}{\sum _{i=-256}^{255}\sum _{j=-256}^{255}{B}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}}\right)}^2+{\left({\displaystyle \frac{\sum _{i=-256}^{255}\sum _{j=-256}^{255}j·B_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}{\sum _{i=-256}^{255}\sum _{j=-256}^{255}{B}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}}\right)}^2}\hfill \\ \hfill & -\sqrt[]{{\left({\displaystyle \frac{\sum _{i=-256}^{255}\sum _{j=-256}^{255}i·G{\left(A\right)}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}{\sum _{i=-256}^{255}\sum _{j=-256}^{255}G{\left(A\right)}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}}\right)}^2+{\left({\displaystyle \frac{\sum _{i=-256}^{255}\sum {j=-256}^{255}j·G{\left(A\right)}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}{\sum _{i=-256}^{255}\sum _{j=-256}^{255}G{\left(A\right)}_{ij}·\mathrm{\Delta }x\mathrm{\Delta }y}}\right)}^2},\hfill \end{array}$$

(10)

where in Eq. (9), i and j denote the sequence numbers of the computing grid in the x and y directions, respectively. $\mathrm{\Delta }x$ and $\mathrm{\Delta }y$ represent the grid sizes in the x and y directions, respectively. The input optical field is denoted as A, and the generator network is represented by G. Meanwhile, B serves as the label corresponding to the intensity distribution generated by numerical simulation. Based on these defined parameters, the modified loss function is formulated as follows:

$$L_G=L_{adv}+\lambda ·L_{pix}+\beta ·L_{BW},$$

(11)

where $L_{adv}$ represents the adversarial loss calculated by the Binary Cross Entropy (BCE) loss function; $L_{pix}$ represents the pixel loss calculated by the $L_1$ loss function. For the discriminator network, we use a BCE loss function corresponding to the adversarial loss. Subsequently, the Adam optimizer is employed with the first-order moment exponential decay rate ${\beta }_1=0.05$ and the second-order moment exponential decay rate ${\beta }_2=0.999$. The hyperparameters $\lambda$ and $\beta$ are set to 5 and 0.5, respectively, based on preliminary experiments to balance the influence of different loss components and optimize the model performance.

3.3. Training Method

To address the GAN’s training challenges, a two-step training method is employed to achieve network convergence gradually. In the first step, a total of 1500 training epochs are conducted. The generator network and the discriminator network are trained alternately at a ratio of 1:1, that is, the generator network is trained once and then the discriminator network is trained once. The labels for the discriminator network are set as follows: the real distribution is labeled as 1, and the fake distribution is labeled as 0. The learning rate of the optimizers for both the generator and discriminator networks is set to 0.001.

In the second step, another 1500 training rounds are carried out. First, the discriminator network is reset to its initial random state. Then, the generator network and the discriminator network are alternately trained at a ratio of 1:2, the discriminator network is trained twice for every one-time training of the generator network. The label setting for the discriminator network is reversed, the real distribution is labeled as 0, and the fake distribution is labeled as 1. The learning rate of the optimizer for the generator network and that for the discriminator network are both set to 0.0001. Throughout the two-step training process, the generator network remains unchanged and is continuously trained. The whole training process utilizes an NVIDIA A100 GPU and takes 4 days to complete.

4. Results

Figure 3 and Figure 4 show the distorted intensity distribution of collimated Gaussian beams propagating through atmospheric turbulence, which are calculated by the generator network and numerical simulation respectively. In these figures, the refractive index structure constant takes values of $1\times {10}^{-16}$, $1\times {10}^{-15}$, and $1\times {10}^{-14}$$m^{-2/3}$, representing weak, moderate, and strong turbulence conditions, respectively. Meanwhile, the transmission distance is set to 2000 m, 4000 m, 6000 m, 8000 m, and 10000 m. From top to bottom, the rows in the figures represent different values of the refractive index structure constant, while from left to right, the columns correspond to different transmission distances.

The figures reveal that the intensity distribution of the laser beam predicted by the neural network aligns well with the outcomes of numerical simulations. Moreover, the radius of the laser beam exhibits an increasing trend with the growth of both transmission distance and turbulence intensity. This observation is in line with the theoretical understanding of laser beam propagation in turbulent media [25]. It strongly suggests that the neural network has effectively captured the influence of turbulence on the propagation of the laser beam.

Similarly, Figure 5 and Figure 6 illustrate the distorted intensity distribution of focused Gaussian beams propagating through atmospheric turbulence, which indicate a high degree of consistency between the predictions of the neural network and the outcomes of numerical simulations. In addition, from the patterns of intensity distribution shown in Figure 6, it is clearly manifested that under weak-turbulence conditions, the focusing spot experiences negligible expansion. However, when a certain turbulence intensity is reached, the converging ability of the Gaussian beam is overshadowed by the expansion effects induced by turbulence, resulting in the formation of a significantly spread spot. This phenomenon is in line with the theory of laser beam propagation in atmospheric turbulence.

By replacing the numerical simulation approach with a neural network, we achieve a significant reduction in computation time, decreasing from 15 seconds to just 6 milliseconds—an impressive reduction of nearly four orders of magnitude. This substantial decrease in computational duration enables the proposed method to effectively meet the demands of applications requiring real-time performance.

To further assess the prediction accuracy of the neural network, we utilized the second moment of the intensity distribution, namely the long-term beam radius, as an index for evaluation. A comparative analysis was conducted between the prediction results obtained from the neural network and those derived from numerical simulations, as illustrated in Figure 7 and Figure 8.

From the results of the long-term beam radius shown in Figure 7 and Figure 8, after being trained on a large-scale training set, the neural network has the ability to simulate the laser beam propagation in the atmospheric turbulent medium, and the relative error under various turbulence intensities and transmission distances does not exceed 16%.

In the case of weak turbulence, both the numerical simulation results and the neural network prediction results exhibit violent oscillations around the theoretical curve, with the overall trend being notably non-smooth. This is because the collimated Gaussian beam in the case of weak turbulence is less affected by turbulence, and its long-term beam radius changes slowly with the increase of the transmission distance. Therefore, even some minor disturbances can cause the statistical results to be magnified in the figure. However, by observing the error magnitude, we can see that although the results seemingly oscillate violently around the theoretical curve, the relative error shows that all the fluctuations are within the error range of 2%, so the actual results are in good agreement with the theoretical curve and numerical simulation results.

Secondly, it is evident that the prediction accuracy of the neural network is higher under weak to moderate turbulence conditions compared to strong and saturated turbulence conditions. This disparity can be attributed to the fragmentation of the laser beam under strong or saturated turbulence, resulting in the formation of numerous small speckles within the two-dimensional intensity pattern. From the perspective of image processing, these speckles manifest as high-frequency signals. The CNN-based neural network encounters challenges in effectively extracting features associated with such high-frequency signals [26], which imposes limitations on its prediction capabilities. Consequently, the error in the neural network’s predictions of intensity distribution is greater in strong and saturated turbulence conditions than in weak and moderate turbulence conditions.

Finally, in the context of saturated turbulence, substantial discrepancies are observed between both the numerical simulation outcomes and the neural network predictions when compared to theoretical expectations. Notably, once the transmission distance exceeds 5000 meters, the theoretical values are significantly elevated relative to those observed in previous cases. This phenomenon primarily arises from the absence of an effective theoretical framework that accurately characterizes the propagation of laser beams under saturated turbulence conditions. The theory employed in this study is primarily designed to address general turbulence scenarios. Thus, the general turbulence theory proves inadequate for accurately modeling such conditions.

5. Conclusion

In this paper, a simulation method for laser beam propagation in turbulent atmosphere based on GAN is proposed. Firstly, a generator based on the U-Net network and a discriminator based on the PatchGAN network are designed. In the construction of the generator, the refractive index structure constant and the transmission distance are introduced as control conditions into the bottleneck layer of the generator. Then, a two-step adversarial training strategy is utilized to complete the training of the generator network, enabling the simulation of Gaussian beams under different turbulence conditions and different transmission distances.

Compared with the traditional numerical method, the neural network method can increase the simulation speed by four orders of magnitude. Furthermore, under various turbulence conditions, the trends of how the second-order statistical moment of the predicted outcomes changes with distance align well with those derived from numerical simulations. Additionally, the relative error is maintained within 16%.

Finally, the differences in the prediction accuracy of the neural network under different turbulence intensities were analyzed. The results show that under weak and moderate turbulence conditions, the neural network prediction results exhibit high accuracy. However, under strong and saturated turbulence conditions, the prediction accuracy of the neural network decreases. This is because the fragmentation of the laser spot leads to the inability of the neural network to effectively capture a large number of high-frequency speckles. Given the complexity of atmospheric turbulence and its potential to introduce intricate patterns, it is important to consider the balance between robustness and accuracy in our models [27,28]. Future research will focus on addressing these challenges to improve prediction accuracy under varying turbulence conditions.