Key Laboratory of Space Laser Communication and Detection Technology, Shanghai Institute of Optics and Fine Mechanics, Chinese Academy of Sciences, Shanghai 201800, China
One fast simulation method using Markov chains was introduced to simulate angular, energy, and temporal characteristics of pulsed laser beam propagation underwater. Angular dispersion of photons with a different number of collisions was calculated based on scattering function and the state transition matrix of Markov chains. Temporal distribution and energy on the receiving plane were obtained, respectively, by use of a novel successive layering model and receiving ratio. The validity of this method was verified by comparing it with the Monte Carlo ray tracing (MCRT) method. The simulation results were close to those obtained by MCRT but were less time consuming and had smoother curves.
The underwater optical wireless communication (UOWC) system makes use of the blue–green wavelength of the visible spectrum and provides higher data rates than the traditional acoustic communication system with significantly lower power consumption over moderate distances[1]. However, the communication performance of the UOWC system depends on the inherent optical properties of water. As water is a kind of multiple scattering medium with intense absorption, when a pulsed laser beam propagates underwater, absorption and scattering will cause angular and temporal spreading, as well as energy attenuation. The temporal spreading, resulting from angular dispersion, limits the reachable data rate, and the energy attenuation limits the range and distance of communication. Thus, for a UOWC system, it is important and necessary to study the performance of a pulsed laser beam through an underwater optical channel. So far, the Monte Carlo ray tracing (MCRT) method[2–8] has been widely used to simulate pulsed laser beam propagation underwater. MCRT is simple to understand and easy to program. Especially, MCRT has been proved to be a reliable method and be in good agreement with experimental results[9–12]. However, the MCRT method needs a long execution time, and the precision is limited by the number of simulated photons.
In this Letter, we firstly introduce Markov chains into the simulation of pulsed laser beam propagation underwater. Using the state transition matrix, angular distribution of photon packets after collision can be easily obtained, as it is only related to the state before the collision[13]. In addition, according to the distribution of the number of collisions and the proposed concept of receiving ratio, we can quickly obtain the energy on the receiving plane. As the state transition matrix can only be used to describe memoryless procedures, but the total propagation delay of one photon packet after multiple scattering is related to the whole propagation path, it is difficult to directly use Markov chains to simulate the temporal distribution of laser pulses on the receiving plane. So, we propose a novel successive layering model to make the Markov chains available for the simulation of temporal distribution of laser pulses. Compared with the MCRT method, the angular distribution after collisions, the energy, and the temporal distribution on the receiving plane simulated by Markov chains have similar results with a higher calculation rate by about two orders of magnitude and a smoother curve. The main reason is that, different from brute-force tracing used in MCRT, Markov chains are based on matrix operation.
The parameters, which relate to pulsed laser beam properties after underwater propagation, mainly include absorption coefficient , scattering coefficient , attenuation coefficient , albedo , and asymmetric parameter . References [14–16] systematically introduced the influences of those parameters and provided a set of data for reference. When a photon undergoes scattering underwater, the Mie theory is available to describe the probability distribution of the scattering angle, provided that the size distribution of scattering particles underwater is known. In simulation, the Henyey–Greenstein (HG) scattering function[17] and Fournier–Forand scattering function[18] have been widely adopted for consistency with experimental data. In this Letter, we select the HG scattering function for its simpler form. The HG scattering function is given by where (ranges from 0 to π) denotes the intersection angle between directions before and after collisions. If setting a local coordinate system that the propagation direction of a photon before scattering is the axis, also denotes the polar of the photon after scattering. The only variable in Eq. (1) is , which is the average value of . Note that, due to symmetry, the azimuthal scattering angle is uniformly distributed over the interval (0–2π). Thus, we just consider the effect of polar angle , and the angular distribution also just refers to .
Sign up for Chinese Optics Letters TOC Get the latest issue of Advanced Photonics delivered right to you!Sign up now
In MCRT, the transmitted laser pulse is treated as a series of photon packets with energy[5]. The interaction of the laser pulse and medium is regarded as collisions between photons and particles, resulting in a change of energy and the propagation direction of the photons. Similar to MCRT, we also treat the absorption and the scattering of laser pulse propagation underwater as a collision procedure between photons and particles. The photons’ energy decay and their directions after collisions follow the HG scattering function. We also use the same terminology adopted in MCRT, such as photons and collisions.
When Markov chains are used to calculate angular distribution, the polar angle should be discredited firstly. We equally divide into parts, with an interval of . In this Letter, the variable is set as 500, since simulation results show that further increasing cannot obviously improve simulation precision but significantly increases computation time. Thus, the angular distribution of a photon can be uniquely denoted by an matrix. According to the HG scattering function, we can form an matrix , i.e., the state transition matrix, to transform the angular distribution before collision to that after collision. The element denotes the probability of a photon with a polar angle ], which scatters to the polar angle after collision. Note that, the state transition matrix refers to the global coordinate system (i.e., the initial laser pulse propagation direction is the axis, and the location laser pulse transmitted into water is original spot), and the HG scattering function refers to the local coordinate system (i.e., the propagation direction of a photon before collision is the axis). We should change the propagation direction of a photon after collision to the global coordinate system when generating the state transition matrix. So, where . Figure 1(a) shows the state transition matrix . It can be seen that the forward scattering is obvious with the elements obviously larger than . Figure 1(b) shows the state transition matrix, raised to the 15th power . It is clear that the angular probability of a photon with the number of collisions of tends to be identical no matter what the initial propagation direction is. The results coincide with Ref. [13].
Figure 1.State transition matrix of water. (a) State transition matrix . (b) State transition matrix raised to the power .
Using the state transition matrix, it is convenient to calculate the angular distribution of a photon with an initial propagation direction along the axis and with the number of collisions , where denotes the initial angular distribution of the photon, which propagates along the axis. Figure 2 shows the angular distribution of a photon with the numbers of collisions of 1, 5, 15, and 30, respectively. The results are compared with those obtained by MCRT. From Fig. 2, it can be concluded that: (1) the results obtained by the state transition matrix are close to those by MCRT; (2) the results obtained by state transition matrix are smooth, but those by MCRT have obvious randomness with many fluctuations (the number of traced photons is one million); (3) as the number of collisions increases, angular distribution tends to the sine function, which stands for totally diffuse light. It should be pointed out that the MCRT is time-consuming, with the computation time about two orders of magnitude longer than that using the state transition matrix.
Figure 2.Normalized angular probability distribution for increasing number of collisions .
Similar to the procedures in MCRT to calculate energy by accumulating the residual energy (due to collision, the photon energy decays) of all photons, which arrive at the receiving plane (due to collision, a portion of photons cannot reach the receiving plane) one by one, we can also obtain the energy on the receiving plane through accumulation by using Markov chains. The difference is that, instead of accumulating energy of the individual photon, we divide the photons into different parts according to and sum up the residual energy of all parts. The procedures are as follows: (1) according to , we divide the photons into a series of parts. The probability of is denoted by the variable , in which denotes the transmission distance. (2) We calculate the receiving ratio , as there are certain probabilities that the photons cannot reach the receiving plane due to collisions. In addition, the energy of a photon decays with a ratio of albedo after collision. (3) Accumulating the residual energy of all the parts, we can obtain the final energy on the receiving plane, with a transmission distance .
To obtain the total energy of all photons arriving at the receiving plane with a certain distance away from the original spot, three parameters are needed; namely, the distribution of the number of collisions , the receiving ratio , and the albedo . In general, is regarded as Poisson distribution with a mean of optical depth [6,7,9], with whereas Eq. (4) is available with the condition of , namely, completely forward scattering. We can still get satisfactory results when forward scattering is intense and optical depth is small. The receiving ratio is defined as the sum of the probability of photons with polar angle , namely
Figure 3 shows the normalized probability distribution of the number of collisions for increasing transmission distance , namely, , where “MCRT” stands for the distributions of photons arriving at the receiving plane obtained by MCRT, “Poisson” stands for Poisson distribution shown in Eq. (4), and “Markov” stands for the product of Markov chains with receiving ratio . It is seen that the distributions of “Markov” are close to those obtained by MCRT. Meanwhile, when transmission increases, the differences within these curves begin to increase.
Figure 3.Normalized probability distribution of the number of collisions for increasing transmission distance .
Finally, the energy of photons arriving at receiving plane is where denotes the energy of the part with , and
We set , , and from Table 1, and in this Letter. The water parameters of Jerlov IB typically correspond to clean water, and those of pool water were gained from actual measurement by the absorption-attenuation spectra (ACS), corresponding roughly to costal water[14,15].
We just consider the photons with the number of collisions less than 30, because the contributions of the photons with the number of collisions larger than 30 can be negligible due to the presence of . Figure 4 shows the energy on the receiving plane versus transmission distance , where denotes the initial energy of the transmitted laser pulse. It is seen that the results obtained by MCRT and Markov chains are roughly identical under both Jerlov IB and pool water, which shows that the Markov chain approach is effective.
Figure 4.Energy loss on the receiving plane versus transmission distance . (a) Jerlov IB: a = 0.06 m−1, b = 0.084 m−1. (b) Pool water: a = 0.117 m−1, b = 0.182 m−1.
Since the angular distribution and the energy obtained by Markov chains are close to those by MCRT, it deserves to be discussed as to whether we can use Markov chains to obtain the temporal distribution of laser pulses on the receiving plane. It is difficult to analyze the temporal distribution directly using the state transition matrix of Markov chains because the angular distribution is only related to the state before collision and has nothing to do with the further earlier states, while temporal distribution needs the total propagation delay from the photon to be transmitted into water until it arrives at the receiving plane. Since the propagation distance of a photon between two consecutive collisions is random, it is hard to calculate the total propagation delay.
To solve this problem, referring to the method used to calculate energy in the previous section, we proposed a novel model, which is named the successive layering model. Next, we introduce, in detail, the procedures to calculate the temporal distribution on the receiving plane by using this model.
The calculation procedures include five steps. (1) Divide the photons into 30 parts according to (). (2) For photons with , divide the transmission distance into layers, which are parallel to each other, following the principle that the collision probabilities of the layers with the same equivalent depth are the same. The equivalent depth is defined as , in which denotes the depth of the layer along with the axis, and denotes the average direction cosine with . , is defined as
Take the situation with the transmission distance and the photons with as an example. Before the first and second collisions, the average direction cosines are and , respectively. In order to satisfy the equation , where and are the depths of the first and the second layer when , and , we set and . In addition, we define the collision spot that is located in the middle of each layer, i.e., the first collision spot is , and the second collision spot is . It should be pointed out that, according to the simulation results, the relative collision location, no matter if it is in the front, middle, or back of each layer, has no obvious influence on the final temporal distribution. Figure 5 shows the transmission procedures of the photon with by using the successive layering model. The photon with is also shown in Fig. 5, in which is the only collision spot. (3) According to the obtained layer depth and angular distribution, calculate transmission distance distribution in the layer. (4) For the photon with , convolute all of the transmission distance distributions obtained in each layer. Subsequently, we can get the temporal distribution of the photons , with and . (5) Sum up all of the distributions of photons with different numbers of collisions according to the energy of each part , and we finally obtain the temporal distribution of all the photons on the receiving plane.
Figure 5.Transmission procedures of the successive layering model, and .
Figure 6 shows the comparisons of the temporal distributions obtained by MCRT and Markov chains, with the transmission distance of 20, 40, 60, and 80 m under Jerlov IB water, respectively. The running time of the MCRT method for 20, 40, 60, and 80 m is about 1562, 3164, 4751, and 6347 s, respectively. Meanwhile, the running time of the less time-consuming Markov chains simulation is 91 s. It is seen that the results obtained by these two methods are very close when the transmission distance is less than 60 m. As is shown in Fig. 6(d), the temporal distribution obtained by Markov chains has a sharper peak and a slower falling edge. The difference between the two curves in Fig. 6(d) mainly comes from the difference of distributions of the number of collisions between Markov chains and MCRT (as shown in Fig. 3). With the increase of transmission distance, the differences of these distributions of the number of collisions increase, which leads to the difference of temporal distributions. The other reason is the simplicity of the successive layering model with very few layers.
Figure 6.Temporal distributions of different transmission distances. (a) , (b) , (c) , (d) .
Though Markov chains have many advantages in the simulation of laser pulses propagation underwater, it still has some drawbacks, such as difficulty to understand and complexity to program. In addition, the energy and spatial distribution in the telescope with a certain field of view and a certain receiving aperture were not investigated yet. Further studies are needed to complete the use of Markov chains in the simulation of laser pulses propagation underwater.
In conclusion, by using Markov chains, we proposed a new method to simulate the characteristics of laser pulse propagation underwater, including angular distribution of photons with a certain number of collisions, the energy, and the temporal distribution on the receiving plane. To simulate temporal distribution, we proposed one novel successive layering model. The results obtained are close to those obtained by MCRT, where they are less time-consuming and have smoother curves. The validity of the method is also verified by experimental results of the pool test. The proposed method combines the advantages of both MCRT and moments technique[19] and provides insights into the field of pulsed laser beam propagation in a multiple scattering medium. Meanwhile, as a newly proposed method, further improvement is needed.