Analysis of Transmission Depth and Photon Number in Monte Carlo Simulation for Underwater Laser Transmission

Abstract

The modeling of laser transmission characteristics in complex seawater is fundamental for bathymetric and oceanographic laser detection systems. Because the factors affecting laser transmission in seawater are independent of one another, firstly, in this paper, a Monte Carlo model of laser propagation in seawaters with suspended matter was established to study the influence of suspended matter with specific radius on the underwater laser transmission. Secondly, the influence of transmission depth and the number of photons on the laser extinction coefficients of seawater containing different concentrations of suspended matter were analyzed, respectively. Thirdly, the relationships between maximum transmission depth, the number of initial photons, and the concentrations were built and verified by simulations. Lastly, an experimental platform was set up and experiments were carried out to verify the Monte Carlo model and the relationships. Results show that (1) both the minimum initial photon number and maximum transmission depth depend exponentially on the concentrations of the suspended matter; (2) the extinction coefficients obtained by the Monte Carlo model and those obtained by experiments are consistent. The absolute values of the differences are less than 0.028 m⁻¹, implying that (1) the proposed Monte Carlo model is effective for simulating laser propagation in seawaters with suspended matter; (2) the established relationships between maximum transmission depth, the minimum initial photon number, and the concentrations of suspended matter have better accuracies, which are valuable for the simulations on attenuation of laser transmission in seawater. The method of this paper can also be extended to the study of suspended solids with other radii and improve the simulation accuracy and decrease simulation time consumption.

1. Introduction

Ocean research is of great significance to Earth science [1]. Generally, marine laser detection systems, such as light detection and ranging (LiDAR), are promising remote sensing technologies for use in oceanographic research, which are widely used in offshore ocean exploration [2], underwater reconnaissance and intelligence gathering [3], mine detection and countermeasure [4], and so on. However, the attenuation of laser in seawater seriously affects the detection performance of marine laser systems. At present, modeling and experimental verification of laser transmission characteristics in complex sea environments are two key technologies in the research of marine laser detection systems. The attenuation consists of Mie scattering from suspended particulate matter, Rayleigh scattering from water and salt molecules, and absorption by dissolved organic matter, phytoplankton, detrital matter, and water. Each part is independent. Generally, attenuation caused by Mie scattering from suspended particulate matter accounts for the largest proportion [5]. Therefore, based on the Mie scattering model, it is particularly important to study the influence of suspended matter on laser transmission in seawater.

Arnush studied the emissivity generated by the initial narrow collimated beam propagating in the ocean based on small-angle scattering theory [6]. Richard et al. adopted Green’s function to solve the radiation transmission equation and obtained the relationship between depth and radiation spatial distribution [7]. Paul et al. established a phenomenological theory of light beam propagation in seawater, and divided the light field into collimated parts and noncollimation parts [8]. Ling et al. used Lambert–Beer’s law and the Mie theoretical model to integrate single matter within the scattering angle range, and obtained the relationship between the large amount of suspended matter contained in oilfield reinjection water and the intensity of the scattered light [9]. The abovementioned methods are numerically accurate but need a huge amount of calculation and are time-consuming; thus, they are difficult to use in some engineering applications that require quick access to computational results.

To obtain results that are consistent with the experimental method and expedite calculation efficiency, semianalytical Monte Carlo methods were used to simulate the optical transmission process in seawater and obtain the extinction coefficient. Zhang et al. combined Mie scattering theory and Monte Carlo numerical simulation to establish an underwater photon transmission model and to analyze the influence of suspended matter on the normalized received energy [5]. Jasman et al. employed the Monte Carlo method to study the scattering characteristics of the underwater diffuse link, and to provide the simulation results of the collimation link for comparison [10]. Liu et al. established a Monte Carlo model for simulating LiDAR signals to evaluate the performance of a shipborne ocean-LiDAR system based on the inherent optical characteristics of seawater and to verify the reliability of the system [11].

However, there are still some problems in Monte Carlo methods for simulating laser transmission in seawater. Li et al. mentioned that each photon in Monte Carlo simulation needs to repeatedly undergo scattering attenuation until it reaches the condition of photon annihilation. Because of the strict limitations of the actual field of view (FOV) and lens aperture, a large number of the photons cannot be received even if a large number of photons are used in the simulation. Therefore, only a very small number of the photons contribute to the received signal and simulation results [12]. As mentioned by Liu et al., a large statistical error inevitably occurs if the number of photons used in simulation is insufficient. Therefore, to decrease the statistical error, it is inevitable to adopt a large number of initial photons for simulation, which requires a relatively long calculation time [11]. Furthermore, Jasman et al. suggested tracking a sufficient number of photons during simulation as much as possible, but did not give a specific initial number of photons [10]. Rafael et al. directly set the number of initial photons as 100 million in their paper [13]. To effectively reduce calculation time on the premise of meeting simulation accuracy, it is necessary to study the most suitable initial photon number for different simulation conditions to reduce the computational complexity and time cost as much as possible without reducing the simulation accuracy.

In this paper, to decide the most suitable initial photon number when simulating laser transmission in seawater containing different concentrations of suspended matter; a Monte Carlo model of laser transmission was established based on Mie scattering theory. The effects of photon number and transmission depth on the extinction coefficient of laser in seawater containing different concentrations of suspended matters were analyzed and verified by experiments.

This paper is organized as follows. The Monte Carlo model is introduced in Section 2. Section 3 analyzes the influence of transmission depth and initial photon number on 532 nm laser extinction coefficients of seawaters containing different concentrations of suspended matter by using the Monte Carlo model. Experimental verification is described in Section 4. Section 5 and Section 6 discuss and summarize the simulation and experimental results.

2. Methodology

In the proposed Monte Carlo modeling, the suspended matter is regarded as spheroidal and the Mie scattering model can be utilized when the radius of the suspended matter is larger than the wavelength of the incident light. The scattering luminous intensity of a suspended particle is denoted as I and described as follows [8,14,15,16,17]:

(1)$I=\frac{{\lambda }^2}{8{\pi }^2{l}^2}{I}_0[i_1(\theta ,\alpha ,n)+i_2(\theta ,\alpha ,n)]$

where $\lambda$ is the wavelength of the incident light; $l$ is the distance between an observation point and the suspended matter; $I_0$ is the luminous intensity of incident light; $i_1$ and $i_2$ are the scattering intensity functions perpendicular and parallel to the scattering surface, respectively, which are determined by scattering angle θ, matter size α, and refractive index n. The two functions are obtained from the Mie scattering theory.

In Monte Carlo modeling, the probability of underwater laser propagation from the transmitter to the receiver is estimated by repeating random photon samplings, which simulate the physical process of laser transmission in seawater containing suspended particles [18,19,20,21]. Firstly, the directions where all photons enter seawater are set, and the lengths of motion of the photons are calculated. Then, the new directions of motion of all photons are calculated by using a scattering phase function. These processes are repeated until all photons annihilate or escape from the boundary. Figure 1 shows the process of Monte Carlo modeling, which consists of four steps and is described as follows:

• ①. Initial setup: The x, y, and z coordinates of the starting point for all photons are set as the origin point (0, 0, 0) of a Cartesian coordinate system, and the positive direction of the z axis is set as the initial transmission direction for all photons. The initial energy of each photon is set as 1.

• ②. Update the position: The transmission step of each photon in seawater is determined by the probability distribution of the photon-free path l, shown as:

  (2)$l=-\frac{\ln r_1}{\sigma }$

  where $r_1$ is the probability of movement of the photon, which is a random number within a uniform distribution from 0 to 1; and $\sigma$ is the extinction coefficient. The azimuth angle of a photon after collision with a particle in seawater is denoted as $\varphi$ and calculated by:

  (3)$\varphi =2\pi r_2$

  where $r_2$ is a random number within a uniform distribution from 0 to 1.

The scattering angle of a photon is denoted as ${\theta }^{*}$ and obtained by numerical statistics based on Equation (4), which satisfies the normalization condition of the scattering angle and best approximates the scattering phase function [14,22].

(4)${\theta }^{\ast }=\{\begin{array}{l}\arccos \left\{\frac{1}{2g}\left[\left(1+g^2\right)-{\left(\frac{1-g^2}{1+g-2g{r}_3}\right)}^2\right]\right\},0\le \theta \le \frac{\pi }{2}\\ \arccos \left\{\frac{1}{2g}\left[{\left(\frac{1-g^2}{1+g-2g{r}_3}\right)}^2-\left(1+g^2\right)\right]\right\},\frac{\pi }{2}\le \theta \le \pi \end{array}$

where g is a parameter used to fit the complex scattering phase function, $\theta$ is the scattering angle of the previous transmission, and $r_3$ is a random number within a uniform distribution from 0 to 1.

The direction cosines of propagation direction of a photon after scattering are denoted as ${\mu }_x^{\ast }$, ${\mu }_y^{\ast }$, and ${\mu }_z^{\ast }$, respectively, and calculated by:

(5)$\{\begin{array}{l}{\mu }_x^{\ast }=\frac{\sin {\theta }^{\ast }}{\sqrt{1-{\mu }_z^2}}({\mu }_x{\mu }_z\cos \varphi -{\mu }_y\sin \varphi )+{\mu }_x\cos {\theta }^{\ast }\\ {\mu }_y^{\ast }=\frac{\sin {\theta }^{\ast }}{\sqrt{1-{\mu }_z^2}}({\mu }_y{\mu }_z\cos \varphi +{\mu }_x\sin \varphi )+{\mu }_y\cos {\theta }^{\ast }\\ {\mu }_z^{\ast }=-\sin {\theta }^{\ast }\cos \varphi \sqrt{1-{\mu }_z^2}+{\mu }_z\cos {\theta }^{\ast }\end{array}$

where ${\mu }_x$, ${\mu }_y$, and ${\mu }_z$ are direction cosines of the photon before scattering. If $\left|{\mu }_z\right|>0.9999$, then ${\mu }_x^{\ast }$, ${\mu }_y^{\ast }$, and ${\mu }_z^{\ast }$ are calculated by:

(6)$\{\begin{array}{l}{\mu }_x^{\ast }=\sin \theta \cos \varphi \\ {\mu }_y^{\ast }=\sin \theta \sin \varphi \\ {\mu }_z^{\ast }=\sin {\mu }_Z\cos \theta \end{array}$

• ③. Update the energy: The energy of each photon decreases after scattering, which is denoted as $W^{\ast }$ and calculated by [23]:

  (7)$W^{\ast }=\frac{\gamma }{\sigma }W$

  where $W$ is the energy of the photon before scattering, and $\gamma$ and $\sigma$ are scattering and extinction coefficients, respectively. In Monte Carlo modeling, it is assumed that the suspended particles are not only equivalent to isotropic and uniform spheres but also evenly distributed in the solution [24]. The scattering coefficient and extinction coefficient are calculated by:

  (8)$\{\begin{array}{l}\gamma =\frac{2N\pi r^2}{{\alpha }^2}{\displaystyle \sum _{k=1}^{\infty }(2k+1)({\left|a_k\right|}^2+{\left|b_k\right|}^2)}\\ \sigma =\frac{2N\pi r^2}{{\alpha }^2}{\displaystyle \sum _{k=1}^{\infty }(2k+1)\mathrm{Re}(a_k+b_k)}\end{array}$

  where N is the number of particles per unit volume, $r$ is matter radius, and $\alpha$ is a dimension parameter and calculated by Equation (9):

  (9)$\alpha =\frac{2\pi r}{\lambda }$

  $a_k$ and $b_k$ in Equation (8) are calculated by:

  (10)$\{\begin{array}{l}{a}_k=\frac{{\psi }_k(a){\psi }_k^{‘}(ma)-m{\psi }_k^{‘}(a){\psi }_k(ma)}{{\xi }_k(a){\psi }_k^{‘}(ma)-m{\xi }_k^{‘}(a){\psi }_k(ma)}\\ b_k=\frac{m{\psi }_k(a){\psi }_k^{‘}(ma)-{\psi }_k^{‘}(a){\psi }_k(ma)}{m{\xi }_k(a){\psi }_k^{‘}(ma)-{\xi }_k^{‘}(a){\psi }_k(ma)}\end{array}$

  where ${\psi }_k\left(x\right)$, ${\xi }_k\left(x\right)$ are calculated by Equations (11) and (12); ${\psi }_k^{‘}\left(x\right)$ and ${\xi }_k^{‘}\left(x\right)$ are derivatives of ${\psi }_k\left(x\right)$ and ${\xi }_k\left(x\right)$, respectively:

  (11)${\psi }_k(x)=\sqrt{\frac{\pi x}{2}}{J}_{k+0.5}(x)$

  (12)${\xi }_k(x)=\sqrt{\frac{\pi x}{2}}{H}_{k+0.5}(x)$

  ${\mathrm{J}}_{k+0.5}\left(x\right)$ and ${\mathrm{H}}_{k+0.5}\left(x\right)$ are the Bessel function of the first kind of half-integer order and the Hankel function of the second kind of half-integer order, respectively [25,26].

  The variable m in Equation (10) is a complex index and calculated by:

  (13)$m=m_1+j{m}_2$

  where $m_1$ and $m_2$ are the real and imaginary parts of the complex refractive index, respectively.

If the energy of a photon is lower than a given threshold, the photon annihilates and is not transmitted. Other photons that have not been annihilated transmit again.

• ④. Photo counting: count the number of photons at a transmission depth $l_z$ along the z direction. Because the photons in any transmission reached the depth $l_z$ should be counted, take the maximum value after comparing with the last statistical value in consideration of the possibility of backscattering. This is to avoid the situation that some photons reach this depth but are not recorded because of backscattering. According to the Beer’s law, the extinction coefficient $\sigma$ is obtained by:

  (14)$\sigma =\frac{d\eta }{d{l}_z},\eta =-\log (\frac{N_r}{N_i})$

  where $N_r$ is the number of the counted photons at $l_z$, and $N_i$ is the number of initial photons.

3. Simulations and Results

3.1. Parameter Setting

• The laser wavelength was set at 532 nm, which is mostly used in LiDAR and the same as it used in experiment. According to the proposed Monte Carlo modeling, the scattering coefficient is related to the complex refractive index of suspended particles in seawater. According to the literature [10], the refractive indexes of the suspended matter are approximately $1.53-\mathrm{j}{10}^{-4}$ and $1.52-\mathrm{j}{10}^{-4}$ at wavelengths of 400 and 600 nm. The complex refractive index was set as $1.52-\mathrm{j}{10}^{-4}$ in the simulation.

• The radius of the active area of the detector used in the simulation was set to 0.1 mm, which is the same as the radius of the commonly used photodetector. For example, the photosensitive surface radius of a widely used photodetector, Thorlabs APD110C, is 100 nm. Because the parameter affected by the receiver radius is the number of received photons, if the minimum number of received photons during simulation can meet the requirements of accuracy, then the simulation results can be obtained with a larger receiver radius.

• The field of view of the detector was set as 90 degrees in the simulation [11,12,13]. According to the definition of scattering extinction coefficient, owing to the scattering of the suspended particles in seawater, part of the radiation deviates from the original propagation direction, resulting in radiation attenuation. The scattering extinction coefficient can be calculated from the relative value of this part of attenuation and transmission depth. All photons scattered in the positive direction to the receiver are regarded as received, and the acceptance range is 180° forward, which is twice the angle of the view parameter.

• Scattering phase function was set to the same function mentioned in [27]. Cornette and Shanks concluded that one-parameter Henyey-Greenstein phase function is nearly consistent with scattering phase function in seawater [28,29,30].

• The radius of the suspended matter in the simulation was set at 3.37 μm, which is the same as $\mathrm{Al}{\left(\mathrm{OH}\right)}_3$ used in experiment, which was provided by Nanjing Senbeijia Biotechnology Limited Company.

• The relationship between the concentration and the number of particles per unit volume in the simulation can be calculated from the Moral formula, which is shown in Table 1.

Based on these parameter settings, the relationship between the extinction coefficient of 532 nm laser and propagation depth were simulated and analyzed.

Figure 2a,b shows the variation curves of the extinction coefficient of 532 nm laser with the transmission depth in seawater containing suspended matter concentrations of 1.5 and 5.0 mg/L, respectively. In the simulation, the number of initial photons is $5\times {10}^5$. The simulated extinction coefficients obtained at different depths fluctuate near a specific value, and the large transmission depth leads to strong fluctuation. Because of the existence of photon annihilation, the number of photons decreases as transmission depth increases. It can be seen that fluctuation of the extinction coefficient rapidly increases when the transmission depth is larger than a threshold. It implies that the number of remaining photons at a large transmission depth is insufficient, and the simulation accuracy is not enough. According to the three-sigma rule, a measurement result is considered exceptional if the fluctuation is three times larger than the standard deviation of all measurement results. As shown in Figure 2a, the corresponding maximum transmission depth is 8 m when the extinction coefficient is within the limited fluctuation. The maximum transmission depth is 1.8 m, as shown in Figure 2b.

Generally, the maximum transmission depth decreases when suspended matter concentration increases. It is necessary to study the relationship between maximum transmission depth and different suspended matter concentrations to obtain a more accurate extinction coefficient by Monte Carlo simulation and to improve simulation speed. In addition, the accuracy of the extinction coefficient obtained by Monte Carlo simulation varies for seawater containing different suspended matter concentrations by using a fixed initial photon number, because the seawater containing a higher suspended matter concentration results in fewer photons remaining to transmit at the same depth. Therefore, the minimum initial photon number for seawaters containing suspended matter concentrations was studied to decrease the simulation error caused by an insufficient initial photon number and to reduce simulation time-consumption caused by an excessive initial photo number. Owing to the randomness of the Monte Carlo method, the extinction coefficients obtained are statistical fluctuations. The extinction coefficient is set as the average of extinction coefficients at different transmission depths, which improves evaluation accuracies of maximum transmission depth and initial photo number.

3.2. Influence of Transmission Depth

The influence of different transmission depths on the extinction coefficients was analyzed. Figure 3a,b shows the variations in extinction coefficients at different transmission depths for seawater containing suspended matter concentrations of 0.9 and 4.0 mg/L, respectively.

It can be seen that the extinction coefficient still fluctuates near a fixed value. For a given concentration, the fluctuation of the extinction coefficient increases with increasing transmission depth. When the transmission depth exceeds a threshold, the extinction coefficient fluctuates strongly, which is because the number of remaining photons is too small. The fluctuation of extinction coefficient is considered reasonable only when the statistical error can be ignored.

Beer’s Law is described by

(15)$P_r=P_z{e}^{-\sigma z}$

where $P_r$ is the received laser power, $P_z$ is the emitted laser power, $\sigma$ is the extinction coefficient, and z is the transmission depth. The emitted power fluctuates up and down within the range of 3% of the accuracy of the commonly used optical power meter Thorlabs PM100D, which will lead to 0.023 ${\mathrm{m}}^{-1}$ of the fluctuation altitude of the extinction coefficient calculated by Equation (15). Furthermore, we can calculate the theoretical attenuation by the empirical formula in [24]. This fluctuation of the extinction coefficient plus the theoretical extinction coefficient is considered as the upper limit, and theoretical extinction coefficient minus the fluctuation of the extinction coefficient is considered as the lower limit. The transmission depth when the fluctuation is beyond the limits is considered as the maximum transmission depth. As shown in Figure 3, the maximum transmission depth in two seawaters is 62 and 5 m. The large differences between the maximum transmission depths for seawater containing different suspended matter concentrations indicate that the same transmission depth cannot be used for studying extinction coefficients in Monte Carlo simulation for seawater containing different suspended matter concentrations; to obtain extinction coefficients with better accuracy, it is important to determine appropriate transmission depths. The maximum transmission depth in 17 seawaters containing different suspended matter concentrations was simulated and calculated, as shown in Table 2.

According to Beer’s law, the number of photons decays exponentially during the transmission process. Therefore, it can be inferred that the maximum transmissible depth also decays exponentially with the increase in suspended matter concentration. To obtain the mathematical relationship between the maximum transmission depth and suspended matter concentrations, five types of exponential functions were used to fit the relationships between the maximum transmission depths and suspended matter concentrations. Table 3 shows the goodness of fit of the five functions.

The third function provides the best goodness of fit, and is shown as:

(16)$y=y_0+A_1(1-e^{-\frac{x}{t_1}})+A_2(1-e^{-\frac{x}{t_2}})$

where $y_0, A_1, A_2, t_1, t_2$ were obtained using least squares. The fitting result follows:

(17)$Z=910.6069-96.01408(1-e^{-\frac{c}{1.40313}})-814.65708(1-e^{-\frac{c}{0.21554}})$

Figure 4 shows the fitting result of the function. We can see that the trend of sampling points and curves is basic () ly the same. These data are the maximum transmission depth at different concentrations when we need to calculate the extinction coefficient by counting the number of received photons at different depths. If the transmission depth is larger than this value, it leads to large fluctuation error and inaccurate results.

The maximum transmission depth can be calculated by Equation (17). It can be seen from Table 4 that the maximum of the relative error is 8%. The relative error is obtained by dividing the difference between the simulation value and the fitting value by the fitting value. It is necessary to take this error into account to calculate the maximum transmission depth. Multiplying the maximum transmission depth in Equation (17) by 0.92 can decrease the fitting error and meet the maximum transmission depth function.

3.3. Influence of Initial Photon Number

In this section, the influence of different initial photon numbers on the mean of extinction coefficients was simulated and analyzed. The extinction coefficients obtained by using different initial photons numbers were analyzed for seawater containing different suspended matter concentrations. Figure 5a,b shows the changes in the simulation extinction coefficients for different initial photon numbers when the concentrations are 0.9 and 4.0 mg/L, respectively.

It can be seen that the fluctuation of extinction coefficients is large for simulations using a small initial photon number. However, initial photon numbers that are too large need larger compute time consumptions and result in low simulation speed. The minimum initial photon number is determined based on the lower and upper limitation. Table 5 shows the minimum initial photon numbers for seawater containing different suspended matter concentrations.

According to Beer’s law, the number of photons in the scattering process decays exponentially. Therefore, five exponential functions were chosen to describe the relationship between the minimum initial photon numbers and suspended matter concentrations. Table 6 shows the fitting results of the five functions.

It can be seen that the goodness of fit of the last function provides the best goodness of fit, but the initial photon numbers are negative at low concentrations. Therefore, the first function is used and shown as:

(18)$y=A(1-e^{-Bx})$

where $A, B$ are obtained by least squares and shown as

(19)$N=-105890(1-e^{0.47813c})$

Figure 6 shows minimum initial photon numbers for seawater containing different suspended matter concentrations. These data are the minimum initial number of photons at different concentrations when we need to calculate the extinction coefficient by counting the number of photons received at different depths. If the initial photon number is less than this value, it leads to large fluctuation error and inaccurate results.

In the simulation of the laser transmission in seawater, the minimum initial photon number can be calculated according to Equation (19). It can be seen from Table 7 that the maximum of the fitting error is 10.8%. The relative error is obtained by dividing the difference between the simulation value and the fitting value by the fitting value. The simulation extinction coefficient can be meet the required fluctuation after calculating the minimum initial photon number multiplied by 1.108.

In reference [24], the laser propagation characteristics in seawater containing suspended solids with radii of 1~4 $\mathsf{\mu }\mathrm{m}$ were simulated, and the number of particles per unit volume is ${10}^9$. According to the Moral law, it can be calculated that the concentration of suspended solids solution is 1.3 × 10⁻¹³ mg/L, and the minimum initial number of photons is 105,890 according to Equation (19). To reduce the influence of fitting error, the simulated initial photon number is set at 117,326 after multiplying by 1.108 for simulation, which is much smaller than the initial photon number of ${10}^6$ given in [24].

Figure 7 shows the normalized received power varies with depths when the initial photon numbers were set at 117,326 and ${10}^6$, as used in [24]. In addition, Table 8 shows the simulation time comparison between the two. Condition 1 is the single simulation time when the initial photon number parameter is set at ${10}^6$, while Condition 2 is the single simulation time when the initial photon number parameter is set at 117,326.

As shown in the legend, four curves show the variation in the normalized received power when the initial photon number is ${10}^6$, as used in [24], and the other four curves show this same function when the initial photon number is 117,326, as calculated by our conclusion. It can be seen that the normalized received power obtained using the two initial photon numbers is consistent, implying that the established relationship is of better accuracy. Furthermore, the simulation with a smaller initial photon number is of faster simulation speed. The simulation time can be reduced to less than 3% by using the conclusion of this study.

4. Experiments and Results

Experiments were carried out to verify the extinction coefficients obtained by the proposed Monte Carlo model and the established relationships [31]. Table 9 shows the simulation extinction coefficients obtained by the proposed Monte Carlo model based on the established relationships between the maximum transmission depth, minimum initial photon numbers, and the suspended matter concentration.

4.1. Experimental Setup

An experiment platform was set up to obtain the extinction coefficient [32,33,34]. Figure 8 shows a schematic diagram of the experimental platform, which consists of a pulsed laser, filter, glass water tank, optical power meter, and data acquisition system. The laser emitted 532 nm laser radiation, which was attenuated by a neutral filter. Then, the attenuated laser passed through a glass tank containing simulated seawater containing different suspended matter concentrations and was detected by an optical power meter. The seawater was simulated using pure water and $\mathrm{Al}{\left(\mathrm{OH}\right)}_3$ suspended matter. The output of the optical power meter was acquired by the data acquisition system. The data obtained were analyzed to calculate the extinction coefficients.

Figure 9 shows the experimental platform [35,36]. Table 10 shows the models of equipment used in the experiment [37,38]. In the experiments, the extinction coefficients of the 532 nm laser for simulated seawater containing different concentrations of suspended matter were measured and compared with simulated results.

4.2. Results

In the experiments, firstly, the power of the laser is 2.3 mW when the glass sink is empty, which is considered as the power of the transmitted laser, excluding the influence of the experimental equipment and air on the extinction coefficient. Secondly, the power of the laser transmitting through the glass tank filled with simulated seawater containing different suspended matter concentrations was measured by the optical power meter and referred as the received powers. The experimental extinction coefficient of the simulated seawater was calculated by Equation (15). Table 11 shows the experimental and simulated and error of extinction coefficients for the simulated seawater containing different suspended matter concentrations.

Figure 10 shows the comparison between the experimental extinction coefficients and those obtained by simulation using the proposed Monte Carlo model based on the established relationships. It can be seen that the two lines are quite consistent. The absolute values of the difference are almost zero, implying that the proposed Monte Carlo model and the established relationship between the maximum transmission depth and the minimum initial photon number were verified by using another method.

5. Discussion

In the previous study on laser transmission simulation, the transmission depth and initial photon number are set by experiences without reference, which affect simulation accuracy and increase simulation time consumption. In this paper, a Monte Carlo model of laser propagation in seawaters with suspended matter was established to study the influence of suspended matter on the underwater laser propagation based on the Mie scattering theory. The following conclusions were reached:

• The transmission processes of 532 nm laser in seawater containing different concentrations of particles were simulated by the Monte Carlo model.

• The effects of different transmission depths and different initial photon numbers on the fluctuation of laser extinction coefficients were analyzed. The upper and lower limits of extinction coefficients were calculated based on the detection accuracy of the commonly used optical powers. The relationship between maximum transmission depths and minimum initial photon number under different concentrations was obtained by a series of Monte Carlo simulations.

• By applying the conclusion of our study to reference [24], it can be found that the initial photon number is greatly reduced, and there is no unacceptable influence on the simulation results. The difference between the two is almost zero. Comparing the simulation time before and after in [24], the simulation efficiency can be improved by 97%.

• Experiments were carried out to compare the extinction coefficients obtained by using the proposed Monte Carlo model and the extinction coefficients obtained by experiments. The two attenuations are basically consistent, implying the two fitting functions described in Equations (17) and (19) obtained in this paper are reliable and can also provide a basis for parameter setting of maximum transmission depth and minimum initial number of photons in future simulations.

6. Conclusions

In this paper, the laser transmission simulation of single-size suspended solids solution is realized, and the simulation and experiment cover the concentration of 0.3–6.5 g/mL. The fitting function can be extended to the study of marine data with lower or higher concentration. The aperture of the receiver set in this simulation is small and the accuracy of the optical power meter is high. The relationship between the extinction coefficient with the maximum transmission depth and the minimum initial photon number can provide guidance for most laser transmission simulations. The theoretical methods in this paper can also be extended to the study of suspended particles with other radii. These works have a certain reference value for the Monte Carlo simulation of laser propagation in the ocean, which can help to improve simulation accuracy and reduce simulation time.

In conclusion, the research of this paper provides two contributions. First, the paper provides a reference to set the transmission distance and the initial number of photons for Monte Carlo simulation of underwater laser transmission. Second, the minimum initial photon number calculated by the fitting function can reduce the simulation time and improve the efficiency while ensuring the accuracy of the result.

Author Contributions

Conceptualization, Y.L. and Y.H.; methodology, Y.L. and D.L.; software, Y.L. and D.L.; validation, Y.L., Y.H. and D.L.; formal analysis, Y.L.; investigation, Y.L.; resources, C.L.; data curation, Y.L.; writing—original draft preparation, Y.L.; writing—review and editing, Y.L.; visualization, D.L.; supervision, Y.H.; project administration, Y.H.; funding acquisition, Y.H. and D.L. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data presented in this study are available on request from the corresponding author.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. The process of Monte Carlo modeling.

Figure 2. Variation curve of extinction coefficients for 532 nm laser with the transmission depths in seawater containing suspended matter concentrations of (a) 1.5 mg/L and (b) 5.0 mg/L.

Figure 3. The mean of the extinction coefficients at different transmission depths for seawater containing suspended matter concentrations of (a) 0.9 mg/L and (b) 4.0 mg/L.

Figure 4. Maximum transmission depth for seawater containing different suspended matter concentrations.

Figure 5. Mean value of the extinction coefficient by using different initial photon numbers when the concentrations are (a) 0.9 mg/L and (b) 4.0 mg/L.

Figure 6. Minimum initial photon number at different concentrations.

Figure 7. Normalized received power varies with depth when the initial photon numbers were set at 117,326 and 1,000,000.

Figure 8. Schematic diagram of the experimental platform.

Figure 9. Experimental platform.

Figure 10. Comparison between experimental extinction coefficient and those obtained by Monte Carlo simulation.

Table 1

The number of particles per unit volume in the simulation for seawater containing different concentrations of the suspended matter used in the experiment.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8
Number of particles per unit volume (×1023)	0.0231	0.0463	0.0695	0.0926	0.1158	0.1389
Concentration (mg/L)	2.1	2.4	2.7	3.0	3.5	4.0
Number of particles per unit volume (×1023)	0.1621	0.1852	0.2084	0.2315	0.2701	0.3087
Concentration (mg/L)	4.5	5.0	5.5	6.0	6.5
Number of particles per unit volume (×1023)	0.3473	0.3859	0.4245	0.4631	0.5017

Table 2

Maximum transmission depths in seawater containing different suspended matter concentrations.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8	2.1	2.4	2.7
Transmission Depth (m)	280	113	62	47	31	26	20	18	15
Concentration (mg/L)	3.0	3.5	4.0	4.5	5.0	5.5	6.0	6.5
Transmission Depth (m)	12	8	5	2	1	3	1	1

Table 3

Goodness of fit for different fitting functions.

Functional Form	$$y=A{e}^{Bx}$$	$$y=e^{A+Bx}$$	$$y=y_0+A_1(1-e^{-\frac{x}{t_1}})+A_2(1-e^{-\frac{x}{t_2}})$$	$$y=y_0+A_1{e}^{-\frac{x}{t_1}}$$	$$y=e^{A+Bx+C{x}^2}$$
Goodness of Fit	0.968	0.968	0.999	0.982	0.978

Table 4

Comparison between the simulation value of maximum transmission depth and the value of the fitting function.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8	2.1	2.4	2.7
Transmission Depth (m)	280	113	62	47	31	26	20	18	15
Fitting Function Value (m)	280	112.9	63	43.9	33.7	26.8	21.5	17.3	14
Relative Error	0	0	0.02	0.07	0.08	0.03	0.07	0.04	0.07
Concentration (mg/L)	3.0	3.5	4.0	4.5	5.0	5.5	6.0	6.5
Transmission Depth (m)	12	8	5	2	1	3	1	1
Fitting Function Value (m)	11.3	7.9	5.4	3.8	2.8	1.9	1	1
Relative Error	0.07	0.02	0.08	0.05	0.05	0.03	0	0

Table 5

Initial photon number for seawater containing different suspended matter concentrations.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8
Initial Photon Number	15,000	30,000	55,000	75,000	105,000	135,000
Concentration (mg/L)	2.1	2.4	2.7	3.0	3.5	4.0
Initial Photon Number	165,000	205,000	255,000	360,000	420,000	600,000
Concentration (mg/L)	4.5	5.0	5.5	6.0	6.5
Initial Photon Number	850,000	1,100,000	1,400,000	1,950,000	2,100,000

Table 6

Average fitting error of different fitting functions.

Functional Form	$$y=A(1-e^{-Bx})$$	$$y=A{e}^{Bx}$$	$$y=e^{A+Bx}$$	$$y=y_0+A_1(1-e^{-\frac{x}{t_1}})+A_2(1-e^{-\frac{x}{t_2}})$$	$$y=y_0+A{e}^{R_{0}x}$$
Goodness of Fit	0.985	0.976	0.976	0.985	0.987

Table 7

Comparison between the simulation value of minimum photon number and the value of the fitting function.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8
Initial Photon Number	15,000	30,000	55,000	75,000	105,000	135,000
Fitting Function Value	16,332	33,183	56,941	82,056	111,044	144,503
Relative Error	0.08	0.1	0.03	0.09	0.05	0.07
Concentration (mg/L)	2.1	2.4	2.7	3.0	3.5	4.0
Initial Photon Number	165,000	205,000	255,000	360,000	420,000	600,000
Fitting Function Value	183,123	227,700	279,151	338,538	458,561	610,996
Relative Error	0.1	0.1	0.09	0.06	0.08	0.02
Concentration (mg/L)	4.5	5.0	5.5	6.0	6.5
Initial Photon Number	850,000	1,100,000	1,400,000	1,950,000	2,100,000
Fitting Function Value	804,599	1,050,490	1,362,780	1,759,400	2,263,140
Relative Error	0.06	0.05	0.03	0.108	0.08

Table 8

The comparison of simulation time.

Radius Size (um)	1	2	3	4
Simulation Time_ Condition 1 (s)	113,450.15	110,598.23	108,943.58	104,684.46
Simulation Time_ Condition 2 (s)	3082.87	2998.25	2625.85	2448.94
Simulation Time_2/Simulation time_1	0.027	0.027	0.024	0.023

Table 9

Simulation extinction coefficient under different concentrations.

Concentration (mg/L)	0.3	0.6	0.9	1.2	1.5	1.8
Simulation Extinction Coefficient (m−1)	0.0511	0.1023	0.1550	0.2047	0.2553	0.3057
Concentration (mg/L)	2.1	2.4	2.7	3.0	3.5	4.0
Simulation Extinction Coefficient (m−1)	0.3581	0.4097	0.4591	0.5118	0.5956	0.6781
Concentration (mg/L)	4.5	5.0	5.5	6.0	6.5
Simulation Extinction Coefficient (m−1)	0.7676	0.8529	0.9347	1.0220	1.1049

Table 10

The models of equipment used in the experiment.

Equipment	Laser	Filter	Optical Power Meter
Model	MPL-H-1064 nm-20 uj-19112652	FL532-3	PD300-TP

Table 11

A side-by-side table with experimental and simulated and error of extinction coefficients of the simulated seawater containing different suspended matter concentrations.

Concentration (mg/L)	0.6	1.2	1.8	2.4	3.0	5.0	6.0
Experimental Extinction Coefficient (m−1)	0.1293	0.2280	0.3115	0.3916	0.4981	0.8809	1.0297
Simulation Extinction Coefficient (m−1)	0.1023	0.2047	0.3057	0.4079	0.5118	0.8529	1.0220
Error Value (m−1)	0.027	0.0233	0.0058	0.0163	0.0137	0.028	0.0077

References

1. Zhang, H.; Zhang, Y. Small Angle Scattering Intensity Measurement by an Improved Ocean Scheimpflug Lidar System. Remote Sens.; 2021; 13, 2390. [DOI: https://dx.doi.org/10.3390/rs13122390]

2. Shimada, S.; Takeyama, Y. Investigation of the Fetch Effect Using Onshore and Offshore Vertical LiDAR Devices. Remote Sens.; 2018; 10, 1408. [DOI: https://dx.doi.org/10.3390/rs10091408]

3. Bikramaditya, D.; Bidyadhar, S. Control of Autonomous Underwater Vehicles: An Overview. Int. J. Autom. Comput.; 2016; 13, pp. 199-225.

4. Bao, N.; Hopkin, D.; Yip, H. Considering Mine Countermeasures Exploratory Operations Conducted By Autonomous Underwater Vehicles. Mil. Oper. Res.; 2014; 19, pp. 19-34.

5. Yingluo, Z. Effect of underwater suspended matters on the transmission characteristics of polarized lasers. J. Opt. Soc. Am.; 2019; 36, pp. 61-70.

6. Arnush, D. Underwater light-beam propagation in the small-angle scattering. J. Opt. Soc. Am.; 1972; 62, pp. 1109-1117. [DOI: https://dx.doi.org/10.1364/JOSA.62.001109]

7. Lutomirski, R.F. An analytic model for optical beam propagation through the marine boundary layer Technical Symposium. Int. Soc. Opt. Photonics; 1978; 160, pp. 110-123.

8. Schippnick, P.F. Phenomenological model of beam spreading in ocean water. SPIE Conf. Ocean. Opt.; 1990; 1302, pp. 13-37.

9. Ling, M.; Fajie, D.; Gencai, S. A concentration measurement model of suspended solids in oilfield reinjection water based on underwater scattering. Measurement; 2018; 117, pp. 125-132.

10. Jasman, F.; Green, R.J. Monte Carlo simulation for underwater optical wireless communications. Proceedings of the 2013 2nd International Workshop on Optical Wireless Communications (IWOW); Newcastle upon Tyne, UK, 21 October 2013; Volume 1, pp. 113-117.

11. Liu, D. Lidar Remote Sensing of Seawater Optical Properties: Experiment and Monte Carlo Simulation. IEEE Trans. Geosci. Remote Sens.; 2019; 57, pp. 9489-9498. [DOI: https://dx.doi.org/10.1109/TGRS.2019.2926891]

12. Jing, L. Monte Carlo study on pulse response of underwater optical channel. Opt. Eng.; 2012; 53, 6001.

13. Mie, G. Beitrage zuroptik trubermedien, speziel kolloidalen metall-losungen. Ann. Phys.; 1998; 25, pp. 377-389.

14. Kraemer, R.M.; Luís, M.; Salgado, H.M. Monte Carlo Radiative Transfer Modeling of Underwater Channel. Wireless Mesh Networks-Security, Architectures and Protocols; IntechOpen: London, UK, 2019; pp. 99-105.

15. Chadi, C.; Mohammj, J. Monte-Carlo-Based Channel Characterization for Underwater Optical Communication Systems. IEEE/OSA J. Opt. Commun. Netw.; 2013; 5, pp. 1-12.

16. Xuewu, D. Propagation and Scattering model of Infrared and Ultraviolet light in Turbid Water. Wirel. Opt. Commun. Conf.; 2013; 6, pp. 601-606.

17. Kaushal, H.; Kaddoum, G. Underwater Optical Wireless Communication. IEEE Access; 2016; 4, pp. 1518-1547. [DOI: https://dx.doi.org/10.1109/ACCESS.2016.2552538]

18. Cox, W.C. Simulation, Modeling, and Design of Underwater Optical Communication Systems; North Carolina State University: Raleigh, NC, USA, 2012; pp. 1-24.

19. Nasir, S. Underwater optical wireless communications, networking, and localization: A survey. Ad Hoc Netw.; 2019; 94, 101935.

20. Guoliang, X.; Da-wei, T. Design of Underwater Laser Triangulation Ranging System and Monte Carlo Simulation Optimization. Adv. Sci. Ind. Res. Cent. Sci. Eng. Res. Cent.; 2018; 7, pp. 1-6.

21. Sanjay, K.S.; Palanisamy, S. A theoretical study on the impact of matter scattering on the channel characteristics of underwater optical communication system. Opt. Commun.; 2018; 408, pp. 3-14.

22. Peppas, K.P.; Boucouvalas, A.C.; Ghassemloy, Z. Performance of underwater optical wireless communication with multi-pulse pulse-position modulation receivers and spatial diversity. IET Optoelectron.; 2017; 11, pp. 180-185. [DOI: https://dx.doi.org/10.1049/iet-opt.2016.0130]

23. Jasman, F. Scattering Regimes for Underwater Optical Wireless Communications using Monte Carlo Simulation. Int. J. Electr. Comput. Eng.; 2018; 8, pp. 2571-2577.

24. Sang, V.Q.; Feng, P.; Tang, B. Study on properties of light scattering based on Mie scattering theory for suspended matters in water. Laser Optoelectron. Prog.; 2015; 52, pp. 231-238.

25. Shen, Z.; Sukhov, S.; Dogariu, A. Monte Carlo method to model optical coherence propagation in random media. J. Opt. Soc. Am. A; 2017; 34, pp. 2189-2193. [DOI: https://dx.doi.org/10.1364/JOSAA.34.002189] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29240093]

26. Kohei, A.; Yasunori, T. Monte Carlo Ray Tracing Simulation of Polarization Characteristics of Sea Water Which Contains Spherical and Non-Spherical Matters of Suspended Solid and Phytoplankton. Int. J. Adv. Comput. Sci. Appl.; 2012; 3, 1313.

27. Shi, Z.; Cui, S.; Xu, Q. Modified analytic expression for the single-scattering phase function. Acta Physica Sinica.; 2017; 68, 180201.

28. Wei, W.; Xiaoji, L. Investigation of 3 dB Optical Intensity Spot Radius of Laser Beam under Scattering Underwater Channel. Sensors; 2020; 20, 422.

29. Yan, D.; Xue, Q. Atmosphere-Ocean Laser Communication Channel Simulation and Modeling. Guangdong Univ. Bus. Stud.; 2009; 4, pp. 554-557.

30. Packard, C.D.; Larsen, M.L. Light scattering in a spatially-correlated matter field: Role of the radial distribution function. J. Quant. Spectrosc. Radiat. Transf.; 2019; 236, 106601. [DOI: https://dx.doi.org/10.1016/j.jqsrt.2019.106601]

31. Yudi, Z. Validation of the Analytical Model of Oceanic Lidar Returns: Comparisons with Monte Carlo Simulations and Experimental Results. Remote Sens.; 2019; 11, 1870.

32. Morel, A. Optical modeling of the upper ocean in relation to its biogenous matter content. J. Geophys. Res. Ocean; 1988; 93, pp. 10749-10768. [DOI: https://dx.doi.org/10.1029/JC093iC09p10749]

33. Roddewig, M.R.; Churnside, J.H.; Shaw, J.A. Mapping the Lidar Attenuation Coefficient in Yellowstone Lake, Yellowstone National Park, USA. Propagation Through and Characterization of Atmospheric and Oceanic Phenomena; Optical Society of America: Washington, DC, USA, 2018; pp. 1-6.

34. Mankovskii, V.I.; Sherstyankin, P.P. Optical characteristics of Lake Baikal waters and their cross-correlations. Izv. Atmos. Ocean. Phys.; 2012; 48, pp. 454-462. [DOI: https://dx.doi.org/10.1134/S0001433812040111]

35. Roddewig, M.; James, H.; Shaw, J.A. Lidar measurements of the diffuse attenuation coefficient in Yellowstone Lake. Appl. Opt.; 2020; 59, pp. 3097-3101. [DOI: https://dx.doi.org/10.1364/AO.389169] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32400591]

36. Hang, L. Iterative retrieval method for ocean attenuation profiles measured by airborne lidar. Appl. Opt.; 2020; 59, pp. C42-C51.

37. André, M.; Stéphane, M. Bio-optical properties of oceanic waters: A reappraisal. J. Geophys. Res.; 2001; 106, pp. 7163-7180.

38. Barth, H. Polychromatic transmissometer for in situ measurements of suspended matters and gelbstoff in water. Appl. Opt.; 1997; 36, pp. 7919-7928. [DOI: https://dx.doi.org/10.1364/AO.36.007919]