We describe the radiative transfer of a Gaussian beam in a water medium using the Monte Carlo method offering basic propagation behaviors. The simulation shows how the energy of the initial Gaussian beam is redistributed as it propagates in coastal water, and also depicts the dependence of the propagation behavior on inherent optical properties of the ocean water such as the single scattering albedo as well as on laser beam parameters, e.g. the M squared. Our results may widen the applicability of LIDARs by providing a couple of design considerations for a bathymetric LIDAR.
I. INTRODUCTION
LIDAR (LIght Detection And Ranging) is a wellestablished technique not only for atmospheric remote sensing and terrestrial hydrography, but also for coastal bathymetry and topography
[1

3]
. In contrast to the wide possible applications, LIDARs are mostly implemented in atmospheric environments in Korea. Only recently great efforts have been paid to developing a laser bathymetry system based on a pulsed green laser with an emphasis on underwater object detection
[4]
.
Numerical techniques for the radiative transfer of an underwater light field have been extensively investigated because understanding the propagation character of a light field in coastal water is a fundamental issue in LIDAR bathymetry as well as in other fields
[5]
. Various models for computing underwater light distributions have been in use, e.g. the discreteordinatemethod
[6]
, the invariant imbedding method
[7]
, and the Monte Carlo method
[8]
. In particular, the Monte Carlo (MC) method has served as a useful tool to model the distribution of underwater downwelling irradiance in the coupled atmosphereocean system, the reflected solar background from the water surface, and the waveform analysis in the timeresolved form
[9]
. In most cases, however, the MC simulation has been used for common ocean optics applications, e.g. the modeling of natural oceanatmosphere environments under a point source like sunlight.
In this paper, we present the radiative transfer of a conventional laser beam, e.g. a Gaussian beam, in a water medium using the Monte Carlo method. The redistribution of a fraction of the photon energy in the original Gaussian beam and in the scattered beam during propagation is described. The dependence of the beam propagation behavior not only on the inherent optical property of the ocean water like the single scattering albedo, but also on the laser beam parameter, e.g. the M squared factor is discussed. In addition, several design considerations for a bathymetric LIDAR are discussed.
II. SIMULATION METHODS
Figure 1
shows the geometrical schematic of the simulation. Photons in a Gaussian distribution start to travel at (
x
,
y
, 0). In a scattering medium with a total attenuation coefficient,
c
(
c
=
a
+
b
, where
a
is the absorption coefficient and
b
is the scattering coefficient), photons are absorbed and scattered as they propagate. At
z
=
z
_{d}
, an arrayed detector is placed to count photons which arrive at the detector grid. The logical flow for calculation is well described in the flow chart in
Fig. 2
. The Mathematica 10.0 program is used for the coding. The spot size of a Gaussian beam at a distance,
z
, is determined by the conventional Gaussian beam theory
[10]
:
Geometrical schematic for the MC simulation.
Flow chart for the MC simulation.
where
w_{0}
is the 1/
e
^{2}
spot size of the laser beam,
z
_{R}
is the Rayleigh length (=
π
w
_{0}
^{2}
/
λ
), and
M
^{2}
is the beam quality factor. The radial position of a photon in a Gaussian distribution is sampled by
[11]
:
where
r
_{i}
and
r
_{f}
are the initial and final radial position, respectively,
q
_{i}
is a random number ranging from 0 to 1. The initial direction cosine is determined by the angle between the
r
_{i}
and
r
_{f}
. For validation, the free space propagation of a Gaussian beam is simulated by MC algorithm and compared to the conventional, showing good agreement.
Once the initial position as well as the direction is determined, a photon travels a geometric pathlength,
s
(=ln(
q
)/
c
), and the position of the photon is updated by:
where (
x
_{0}
,
y
_{0}
,
z
_{0}
) is the previous coordinate, (
μ
_{x}
,
μ
_{y}
,
μ
_{z}
) is the direction cosine. After a scattering event, the direction cosine should be appropriately updated. For validation of the direction update, the diffusion of a point source in an isotropically scattering medium is simulated yielding a good agreement with the diffusion theory. The scattering character in this study is determined by the HenyeyGreenstein phase function
[12]
:
where
g
is the asymmetry factor. For coastal water, the asymmetry factor is taken to be 0.9
[13]
. In each scattering event, the photon weight should also be renewed by:
where
ω
_{0}
=(
c
−
a
)/
c
is the single scattering albedo. The parameters in use for the simulation are listed in
Table 1
.
The parameters in use for the MC simulation
The parameters in use for the MC simulation
III. RESULTS AND DISCUSSION
The inset of the
Fig. 3(a)
shows a density plot of the laser beam and its projected profile at z = 1.1
z
_{R}
(~ 26 m). There is a strong peak in the center, which is from the original Gaussian beam, while the scattered photons are distributed around the center peak. As the laser beam propagates in a water medium, the photons are scattered and the width of the scattered photon distribution is widened as shown in
Fig. 3(a)
. The 1/
e
^{2}
width of the scattered photon distribution (red dots) is deduced by a Gaussian fitting of the projected beam profile. For the given condition in the
Table 1
, the width grows as large as 10.8 m at
z
=
z
_{R}
(~ 23.6 m), which far exceeds the original beam width. The widths deduced by the MC method are compared with the results of Lutomirski’s theory (blue line)
[14
,
15]
, yielding good agreement. For calculation by the analytic model, the waveslope variation (
σ
_{x,y}
), rms scattering angle owing to maritime aerosols (
σ
_{β}
) are neglected and the incidence angle is set to 0, and the rootmeansquare scattering angle (
Θ
_{R}
) is calculated by the asymmetry factor
g
. The diffuse attenuation coefficient (
K
_{d}
) is deduced by inherent optical properties in
Table 1
[16]
.
Figure 3(b)
shows the variation of the fraction of photons in the center peak (red dots) and in the scattered distribution (blue rectangles) at various depths. As the laser beam propagates, the fraction of photons in the center peak decreases, while the energy in the scattered distribution increases. This can be understood in terms of the redistribution of the energy in the center peak into that in the scattered distribution through the diffusion. After the beam reaches a depth of
z
= 1.6
z
_{R}
(~ 38 m), the center peak completely disappears, thus only the fraction of the photons in the scattered distribution propagates. The conventional amplitude of the return signal from a target, e.g. sea bed and an underwater object is mainly contributed by the original Gaussian photon fraction because the scattered photons are scattered out. Thus the maximum detectable depth using a LIDAR return signal can be expected to be more or less 38 m in the given environmental condition. On the other hand, the maximum detectable depth can be further extended by increasing the output energy of the laser and detecting tiny a fraction of the photons in the scattered distribution using a singlephoton detector, e.g. a photomultiplier tube. In this case, the fieldofview loss factor induced by the broad distribution of scattered photons should be appropriately estimated for the optimized LIDAR performance
[15
,
17]
.
(a) The 1/e^{2} width of the scattered fraction as a function of the depth deduced by MC method (red dots) and Lutomirski model (blue line). Inset, a density plot of the laser beam and its projected profile at z = 1.1z_{R} (~ 26 m). (b) The fraction of photons in the center peak (red dots) and in the scattered distribution (blue rectangles) as a function of the depth.
The width of the fraction of the photons in the scattered distribution as a function of the single scattering albedo (
ω
_{0}
) at
z
=
z
_{R}
is shown in
Fig. 4(a)
. For the calculation, the water absorption is fixed to be 0.05/m. The influence of the single scattering albedo on the beam width of the scattered fraction is not critical. The standard deviation of the various beam widths is only 6.5 % of the mean width. Thus, the environmental parameter is not the major factor in the beam spread of the laser beam in the scattering medium.
Figure 4(b)
shows the energy redistribution behavior as a function of the single scattering albedo. The fraction of photons in the original Gaussian beam (red dot) starts to significantly decrease at
ω
_{0}
~ 0.33, while the fraction of photons in the scattered distribution (blue rectangles) starts to increase at this point. It can be easily expected that the return signal amplitude starts to be degraded at
ω
_{0}
~ 0.33 in the given laser output energy.
(a) The 1/e^{2} width of the scattered fraction, and (b) the fraction of photons in the center peak (red dots) and in the scattered distribution (blue rectangles) as a function of the single scattering albedo.
In
Fig. 5(a)
, the dependence of the scattered beam width on the beam quality factor for
z
=
z
_{R}
is shown. For various beam quality factors ranging from the ideal case (1.0) to the multimode case (20.0), the scattered beam width is almost constant at ~ 10.9 m. Furthermore, the fraction of photons, i.e. the averaged energy density has no influence on the beam quality factor as shown in
Fig. 5(b)
. The photon fraction in the Gaussian beam and the scattered distribution are almost constant at 0.13 and 0.87, respectively. This result shows that the propagation behavior in terms of the averaged energy density in the original Gaussian beam and the scattered photon distribution is independent of the Msquared. This is true even if the fine structure of the Gaussian beam is distorted when taking refractive index fluctuations of the coastal water as well as the random evolution of the ocean surface into account. These factors are not considered in the current calculation
[13
,
17]
. Thus any laser source with high beam quality is not required because the initial beam quality diminishes during the propagation after a distance, e.g. longer than 38 m in the given environmental condition. In a practical point of view, the effort required to maintain single mode character in high power laser for a bathymetric LIDAR system can be significantly reduced.
(a) The 1/e^{2} width of the scattered fraction, and (b) the fraction of photons in the center peak (red dots) and in the scattered distribution (blue rectangles) as a function of the beam quality factor, the M squared.
IV. CONCLUSION AND OUTLOOK
In conclusion, the Monte Carlo simulation for underwater propagation of a Gaussian beam was performed and the results were presented. The variation of a fraction of the photon energy in the original Gaussian beam and in the scattered beam during propagation was simulated. The dependence of the beam propagation behavior on the inherent optical properties of the ocean water like the single scattering albedo and on the laser beam parameters, e.g. M
^{2}
factor was discussed. The results may offer several design considerations of a bathymetric LIDAR, e.g. the fieldofview loss factor, the requirement of single mode character of the laser.
The simulation can be further improved by including the random ocean surface, the measured volume scattering function, and the measured IOPs, e.g. the absorption coefficient, the single scattering albedo.
Guenther G. C.
,
Thomas R. W. L.
,
LaRocque P. E.
1996
“Design considerations for achieving high accuracy with the shoals bathymetric lidar system,”
Proc. SPIE
2964
26 
37
Guenther G. C.
,
Cunningham A. G.
,
LaRocque P. E.
,
Reid D. J.
2001
“Meeting the accuracy challenge in airborne LIDAR bathymetry,”
Proc. EARSeLSIGWorkshop LIDAR
Dresden, Germany
1
1 
27
Wang C. K.
,
Philpot W. D.
2007
“Using airborne bathymetric lidar to detect bottom type variation in shallow waters,”
Remote Sens. Environ.
106
123 
135
DOI : 10.1016/j.rse.2006.08.003
Kim J. I.
,
Yoo H. K.
,
Cho J. J.
,
Kim H. R.
,
Kim W. H.
,
Jeon Y. G.
2014
“A green laser system for detecting underwater objects and its applications,”
Proc. The 10th Marine Weapon Conference
Agency for Defense Development, Jinhae, Korea
101 
102
Jerlov N. G.
1976
Marine Optics
Elsevier
Amsterdam, Netherlands
Stamnes K.
,
Dale H.
1988
“Numerically stable algorithm for discreteordinatemethod radiative transfer in multiple scattering and emitting layered media,”
Appl. Opt.
27
2502 
2509
DOI : 10.1364/AO.27.002502
Mobley C.
1989
“A numerical model for the computation of radiance distributions in natural waters with windroughened surfaces,”
Limnol. Oceanogr.
34
1473 
1483
DOI : 10.4319/lo.1989.34.8.1473
Gordon H.
,
Brown O.
,
Jacobs M.
1975
“Computed relationships between the inherent and apparent optical properties of a flat homogeneous ocean,”
Appl. Opt.
14
417 
427
DOI : 10.1364/AO.14.000417
Guenther G. C.
,
Maune D. F.
2001
Digital Elevation Model Technologies and Applications: The DEM User’s Manual
ASPRS
“Airborne lidar bathymetry,”
Siegman A. E.
1986
Lasers
University Science Books
Milsom P. K.
2000
“A rayoptic Monte Carlo, description of a Gaussian beam waistapplied to reverse saturable absorption,”
Appl. Phys. B
70
593 
599
DOI : 10.1007/s003400050867
Henyey L. C.
,
Greenstein J. L.
1941
“Diffuse radiation in the galaxy,”
Astrophys. J.
93
70 
83
DOI : 10.1086/144246
Kirk J. T. O.
1981
“Monte Carlo procedure for simulating the penetration of light into natural waters,”
Theor. Appl. Genet.
60
197 
214
DOI : 10.1007/BF02342540
Lutomirski R. F.
1978
“An analytical model for optical beam propagation through the marine boundary layer,”
Proc. SPIE
160
110 
120
Tulldahl H. M.
,
Steinvall K. O.
1999
“Analytical waveform generation from small objects in lidar bathymetry,”
Appl. Opt.
38
1021 
1039
DOI : 10.1364/AO.38.001021
Guenther G. C.
1985
“Airborne laser hydrography: system design and performance factors,” NOAA professional paper series
National Ocean Service
1 
385
Carr D.
,
Tuell G.
2014
“Estimating fieldofview loss in bathymetric lidar: application to largescale simulations,”
Appl. Opt.
53
4716 
4721
DOI : 10.1364/AO.53.004716
Mahdieh M. H.
2008
“Numerical approach to laser beam propagation through turbulent atmosphere and evaluation of beam quality factor,”
Opt. Commun.
281
3395 
3402
DOI : 10.1016/j.optcom.2008.02.040
Smith R. C.
,
Baker K. S.
1981
“Optical properties of the clearest natural waters (200800 nm),”
Appl. Opt.
20
177 
DOI : 10.1364/AO.20.000177