© The Author(s), published by EDP Sciences, 2026

1 Introduction

Ocean plays a vital role in regulating the global climate, sustaining biodiversity, and supporting essential resource [1–3]. Ocean sensing is therefore a fundamental prerequisite for understanding and sustainably utilizing the ocean. However, conventional detection methods face intrinsic limitations: sonar systems aren’t able to detect through air–water interface, while passive satellite remote sensing is not able to provide the vertical information on the upper ocean layers [4–7]. Among ocean sensing techniques, Light Detection and Ranging (Lidar) has emerged as an effective alternative. As an active remote sensing system, lidar penetrates the air–water interface and captures high spatial and temporal resolution profiles of subsurface structures [8–10]. These capabilities establish lidar as an essential complement to existing ocean sensing approaches.

Lidars has been developed for diverse applications, among which one of the earliest and most prominent uses is topobathymetry [11]. Modern topobathymetry lidar systems are capable of achieving high-precision integrated coastal and terrestrial topographic mapping [12–15], seabed substrate classification [10, 16], as well as surveying in turbid waters and shallow areas [17]. Moreover, lidar profiling information provides essential support for studying marine organisms and retrieving seawater optical properties. Investigating the distribution and dynamics of phytoplankton [18–22] based on lidar profiling data contributes to a deeper understanding of the mechanisms driving marine primary productivity and the global carbon cycle. Lidar has also been proven to be highly effective in surveying fishery resources [23–27] near the sea surface, and its performance on deriving optical properties [28–32], such as the attenuation coefficient and volume backscattering, has been validated by shipborne in-situ measurements. Additionally, lidar enables measurements of wave characteristics without contacting the water body or disturbing hydrodynamic processes, serving as an effective tool for observing ocean surface waves [33, 34] and internal waves [35, 36].

Underwater target detection is now emerging as a new direction for lidar applications, which has importance influences on human activities, such as navigation planning and port construction [2, 23]. In the field of underwater target detection, current research mainly focuses on two directions: underwater imaging and rapid detection of underwater targets. Since the first application of single-photon avalanche diode (SPAD) arrays to underwater imaging is represented, single-photon imaging technology has steadily advanced [37]. Subsequently, the first fully submerged underwater lidar transceiver system based on single-photon detection technologies has been demonstrated, achieving images of stationary targets at ranges up to 7.5 attenuation lengths [38]. More recently, a compact, all-fiber underwater imaging lidar has been presented that achieved imaging distances exceeding 10 times the optical attenuation length, further highlighting the progress in underwater imaging lidar research [39]. However, the application of lidar for rapidly detecting underwater targets from above water surface remains limited. Only some mapping lidar systems possess the capability for rapid target detection, and there is no dedicated lidar system specifically designed for this purpose. The underwater obstruction detection capabilities of the Scanning Hydrographic Operational Airborne Lidar Survey (SHOALS) were discussed by adopting the enhanced object detection algorithm [12, 13, 40]. A lightweight UAV-borne lidar Mapper4000 U was developed for shallow water mapping, which had the ability to recognize the existence of a 1-m cube and the rough shape of a 2-m cube at a depth of 12 m by point cloud fitting [41].

There is also a great deal of algorithms to generate lidar signals, such as the lidar equation, analytical model and Monte Carlo simulation technique. The lidar equation is able to describe depth-dependent lidar signals quickly but has limitations in dealing with multiple scattering [42, 43]. The analytical model is highly difficult to solve due to its complexity. Only approximate solutions can be achieved through some idealized assumption [44]. The standard Monte Carlo method is based on a huge number of random calculations of photons’ movements in the water, and photons can be accepted only when they finally arrive at the receiver aperture. The condition is too strict to reconstruct the lidar profiles, consuming a large amount of computing power and time [45]. The semi-analytic Monte Carlo radiative transfer model for oceanographic lidar systems (SALMON) is based on the method of stochastic and analytic techniques, which estimates the contribution of each scattered photon by the function from a particle scattering position to the detector [46]. With the advantages in computational efficiency and accuracy, the SALMON model has been the foundation of the MC simulator for oceanic lidar. After decades of development, the semi-analytical MC model has been applied in seabed mapping, retrievals of bio-optical properties and polarized signal analysis and so on [47–52]. However, current research in ocean mapping or target detection is primarily focused on the seabed modeling and has some overly idealized assumptions, which assumes that the seabed fills the entire lidar field of view and simplifies the photon–seabed interactions, without accounting for three-dimensional interactions or the blocking effect on the beam transmission [53, 54]. Besides, the majority of these studies assume a homogeneous water body and neglect environmental noises, making it difficult to capture the complexity of realistic detection scenarios. On the other hand, the simulators for oceanic lidar profiling have made some advancements in accuracy. While these simulators mainly focus on the slope variations of lidar returns that are closely related to signal retrieval, photon propagation only needs to account for the effects of absorption and scattering by the water, which makes them incapable of detecting three-dimensional target.

As it stands, lidar systems that are applicable to or specifically designed for underwater target detection have yet to be widely developed. A dedicated Monte Carlo simulator for three-dimensional underwater target detection has yet to be established. Three-dimensional simulations are crucial for capturing spatial heterogeneity, beam divergence, and realistic scanning geometries. Therefore, it is highly necessary to develop a simulator that not only emphasizes accuracy but also generates three-dimensional target detection as well as a lidar system designed base on the simulator’s guidance.

Our research group has already conducted extensive prior work in this area. A standard MC model to detect infinite planar target [55] and a semi-analytical MC model to detect a two-dimensional flat plane with limited size [56] have been proposed. In this paper, we present a refined semi-analytical Monte Carlo model for underwater three-dimensional target detection (MCT), especially including the interaction process between the beam and the target, the blocking effect on the beam transmission and the reception of the target reflection signals. The MCT model can simulate the scanning signals for underwater three-dimensional targets through surface waves in stratified water columns, and incorporate realistic environmental noise to better approximate field conditions. We develop a miniaturized, linear scanning oceanic lidar system (SOL) deployable on commercial UAV platforms. The MCT model is validated by comparisons with the lidar equation results and field measurements from SOL. Furthermore, using MCT, we evaluate SOL’s detection capabilities in Jerlov II and Jerlov 3C waters, define an extended detection range metric, and examine its relationship with vertical penetration and horizontal scanning resolution, offering an integrated framework for system performance assessment and design optimization in realistic marine environments.

In the rest of this paper, we present the MCT model in Section 2. We present the SOL system in Section 3. The MCT model is validated in Section 4. The target detection capabilities of SOL in different waters are evaluated in Section 5. The conclusion is made up in Section 6.

2 Monte Carlo model of underwater target detection

For a complete and detailed introduction to the MCT model, a cylinder target is chosen as a representative example. The cylindrical geometry captures essential features relevant to underwater lidar detection, including extended surfaces, well-defined edges, and stable occlusion characteristics under scanning illumination. As a representative intermediate case, the cylindrical geometry combines manageable geometric complexity with physically meaningful target characteristics relevant to practical detection scenarios. It is worth emphasizing that the proposed MCT model is formulated in a geometry-independent manner. Other target can be incorporated by modifying the geometric description.

The main parts are introduced below, and the flow chart of the model is shown in Figure 1. The whole simulation starts with the emission of photons and ends with their final death. During the propagation process, the photon can be scattered, reflected by the target, escaping from the boundary, or received by the detector. The surface of the target is seen as ideal Lambertian. The inherent optical properties (IOPs) model of water used in simulation is based on previous research [43, 57–62].

Fig. 1 Flow chart of MCT simulation model.

2.1 Launching

The Cartesian coordinates of the simulation model are shown in Figure 2a, the initial position of the lidar system is P_st  = [x_st, y_st, 0], and the downward direction of the z-axis is positive. The intensity of the emitted beam follows a Gaussian distribution, and the initial position of the photon is sampled from the laser spot on the water surface (Fig. 2b). The intensity distribution on the water surface is [63]$$\begin{array}{c}I\left(\mathrm{ }{R}_{\mathrm{G}}\right)=\frac{I_0}{\sqrt{2\pi }\sigma }\exp \left(-\frac{R_G^2}{2{\sigma }^2}\right),\end{array}$$(1)

Fig. 2 (a) Global coordinate system of the simulation model, and (b) sampling the point from the Gaussian spot.

where σ is standard deviation of Gaussian distribution. The sampled spot radius R_G is expressed as$$\begin{array}{c}{R}_G=\sigma \sqrt{-2\ln \left({\xi }_{\mathrm{rnd}}\right)},\end{array}$$(2)

where ξ_rnd is a random number uniformly distributed over (0, 1). The initial azimuth angle on the Gaussian spot is sampled as follow$$\begin{array}{c}{\phi }_0=2\pi \bullet {\xi }_{\mathrm{rnd}}.\end{array}$$(3)

The photon’s initial position on the water surface then be obtained as$$\begin{array}{c}\{\begin{array}{c}x=R_G\cos {\phi }_0+x_{\mathrm{st}}\\ y=R_G\sin {\phi }_0+y_{\mathrm{st}}\\ z=H \end{array},\end{array}$$(4)

where H is the flight height. The lidar system positions x_st and y_st will be changed with the UAV moves. The initial incident direction of the photon is$$\begin{array}{c}\{\begin{array}{c}{U}_0\left(x\right)=\sin {\theta }_0\cos {\phi }_0\\ U_0\left(y\right)=\sin {\theta }_0\sin {\phi }_0\\ U_0\left(z\right)=\cos {\theta }_0\end{array},\end{array}$$(5)

where the emitting angle can be calculated by θ₀ ≈ tan⁻¹(R_G/H). When the photon passes through the air–sea interface, the instantaneous state of the waves must be taken into account. The normal vector of the surface wave is defined as$$\begin{array}{c}\{\begin{array}{c}{U}_w\left(x\right)=\sin {\theta }_w\cos {\phi }_w\\ U_w\left(y\right)=\sin {\theta }_w\sin {\phi }_w\\ U_w\left(z\right)=\cos {\theta }_w\end{array},\end{array}$$(6)

where ${\phi }_w$ is also sampled by equation (3). The normalized probability density function of θ_w is proposed by Cox and Munk [64]:$$\begin{array}{c}p\left({\theta }_w\right)=\frac{2}{{\phi }_w^2}\exp \left[\frac{-{\tan }^2\left({\theta }_w\right)}{{\sigma }^2}\tan \left({\theta }_w\right){\sec }^2\left({\theta }_w\right)\right],\end{array}$$(7)

where σ is the root-mean-square of surface slope and calculated by $\sigma ={\left(0.003+0.00512\nu \right)}^{1/2}$. ν is the wind speed. The propagation direction U of the photon after entering the water can be obtained from U₀ and U_w by applying Snell’s law.

2.2 Moving and scattering

Absorption and scattering events happen when the photon propagates underwater. According to the randomly sampled optical distance, the photon’s moving direction, and the attenuation coefficient profile of water, the geometric distance of the photon’s random movement can be calculated [65]. The geometric distance is used to update the photon’s position. At the new position, the photon’s weight will be reduced, and the moving direction will be changed. For the scattering direction, its scattering polar angle θ can be sampled by the scattering phase function, while the scattering azimuth angle φ is uniformly distributed over [0, 2π]. The above process can be found in the reference [49]. However, considering that the concentration of scattered particles in the stratified water varies with depth, which will affect the probability of the photon’s colliding with particles, this effect must be included in the simulation. In this study, the scattering phase function of water molecules mixed with particles is used for scattering polar angle θ sampling$$\begin{array}{c}\tilde{\beta }\left(\theta ,z\right)=\frac{b_w}{b}{\tilde{\beta }}_w\left(\theta \right)+\frac{b_p}{b}{\tilde{\beta }}_p\left(\theta ,z\right),\end{array}$$(8)

where the ${\tilde{\beta }}_w\left(\theta \right)$ is the normalized phase function of water and ${\tilde{\beta }}_p\left(\theta \right)$ is the Petzold average-particle phase function [66–68]. b is the total scattering coefficient, which is the sum of water scattering coefficient b_w and the particles scattering coefficient b_p respectively. It should be noted that, the scattering polar angle θ is sampled with respect to the current photon propagation direction U, which serves as the reference axis of the local coordinate system. Consequently, to obtain the updated photon propagation direction U_new, the directional transformation from the local scattering coordinate system to the global coordinate system is given by$$\begin{array}{c}\left[\begin{array}{c}{U}_{\mathrm{new}}\left(x\right)\\ U_{\mathrm{new}}\left(y\right)\\ U_{\mathrm{new}}\left(z\right)\end{array}\right]=\left[\begin{array}{ccc}\frac{U\left(x\right)U\left(z\right)}{\sqrt{1-{U\left(z\right)}^2}}& \frac{-U\left(y\right)}{\sqrt{1-{U\left(z\right)}^2}}& U\left(x\right)\\ \frac{U\left(y\right)U\left(z\right)}{\sqrt{1-{U\left(z\right)}^2}}& \frac{U\left(x\right)}{\sqrt{1-{U\left(z\right)}^2}}& U\left(y\right)\\ -\sqrt{1-{U\left(z\right)}^2}& 0& U\left(z\right)\end{array}\right]\left[\begin{array}{c}\sin \theta \cos \phi \\ \sin \theta \sin \phi \\ \cos \theta \end{array}\right].\end{array}$$(9)

2.3 Interaction with the target

The target underwater is set up as a cylinder with a longitudinal length L_C and a cross-sectional radius R_C. As shown in Figure 3, the photon’s position is P_i, based on the sampled geometric distance L_step and the incident direction U, and the new position P_step can be calculated by P_step = P_i + L_step· U. If the position P_step(z) ≥ H + d_w − R_C, and the moving direction U(z) > 0, this will indicate that the photon moves downwards and it may collide with the target. When U(z) < 0 and P_step(z) ≤ H + d_w + R_C, it indicates a possible upward interaction event. Assume L_touch is the distance between P_i and P_touch, and the interaction position P_touch can be described by P_touch = P_i + L_touch· U. According to the projection of the ray interaction with the cylinder in Figure 3, P’_C is the center of the projection circle, and L’_foot is the distance from P’_i to P’_foot. Since P’_i = [y_i, z_i] is known in the yoz plane, the point P’_foot can be represented as [y_i + L’_foot·U(y), z_i +L’_foot·U(z)]. Based on vector $\stackrel{⃗}{{\mathit{P}}_{\mathit{C}}\mathit{\text{'}}{\mathit{P\text{'}}}_{\mathbf{foot}}}$ is perpendicular to vector $\stackrel{⃗}{{\mathit{P}}_{\mathit{i}}\mathit{\text{'}}{\mathit{P\text{'}}}_{\mathbf{foot}}}$, the distance L’_foot and the point P’_foot can be obtained by$$\begin{array}{c}\stackrel{⃗}{{\mathit{P}}_{\mathit{C}}\text{'}{{\mathit{P}}^{\text{'}}}_{\mathbf{foot}}}\bullet {\vec{\mathit{e}}}_{\mathit{i}}=0,\end{array}$$(10)

Fig. 3 Geometry of the photon colliding with the cylinder target, the left sub-graph is projection of the interaction process.

where ${\vec{\mathit{e}}}_{\mathit{i}}=\left[\mathit{U}\right(y),\mathit{U}(z\left)\right]$ is the unit vector from $\stackrel{⃗}{{\mathit{P}}_{\mathit{i}}\mathit{\text{'}}{\mathit{P\text{'}}}_{\mathbf{foot}}}$. For an interaction event, the distance $\left|\stackrel{⃗}{{\mathit{P}}_{\mathit{C}}\mathit{\text{'}}{\mathit{P\text{'}}}_{\mathbf{foot}}}\right|$≤ R_C (radius of the cross section). The distance $L_{\mathrm{touch}}^{\text{'}}$ is expressed as$$\begin{array}{c}{L}_{\mathrm{touch}}^{\text{'}}=L_{\mathrm{foot}}^{\text{'}}-\sqrt{R_c^2-{\left|\stackrel{⃗}{{\mathit{P}}_{\mathit{C}}\mathit{\text{'}}{{\mathit{P}}^{\mathit{\text{'}}}}_{\mathbf{foot}}}\right|}^2}.\end{array}$$(11)

$L_{\mathrm{touch}}^{\text{'}}$ in projecting space is converted to three-dimensional space by$$\begin{array}{c}{L}_{\mathrm{touch}}=\frac{L_{\mathrm{touch}}^{\text{'}}}{\sqrt{{\mathit{U}\left(y\right)}^2+{\mathit{U}\left(z\right)}^2}},\end{array}$$(12)

and the coordinates of P_touch then can also be calculated. If P_touch(x) $\in \left[\frac{-L_C}{2},\frac{L_C}{2}\right]$ and the geometric distance of this movement isn’t less than L_touch, the collision eventually occurs.

To calculate the direction after collision, the differential tangent plane around ${\mathit{P}}_{\mathbf{touch}}^{\mathit{\text{'}}}$ is treated as a Lambertian reflector, so the probability density function of the reflected angle θ_r can be expressed as$$\begin{array}{c}P\left({\theta }_r\right)=\cos {\theta }_r,\end{array}$$(13)

then the propagating direction Ut after collision can be obtained by substituting θ_r into θ of (7) and the unit direction vector of $\stackrel{⃗}{{\mathit{P}}_{\mathit{C}}\mathit{\text{'}}{\mathit{P\text{'}}}_{\mathbf{touch}}}$ into U, where the reflected azimuth angle φ_r = 2πξ_rnd. For the upward interaction event, U_t(z) takes a positive value; for the downward interaction event, U_t(z) takes a negative value. The photon’s weight is reduced to W_i = ρW_i, and ρ is the surface reflectivity of the target.

2.4 Interaction with the target

Only if the photon’s position P_i is in the volume of detector’s field of view (FOV), the signal contribution can be recorded. The condition of the detection is expressed as [69]$$\begin{array}{c}\{\begin{array}{c}\sqrt{{\left[{\mathit{P}}_{\mathit{i}}\left(x\right)-{\mathit{P}}_{\mathbf{FOV}}\left(x\right)\right]}^2+{\left[{\mathit{P}}_{\mathit{i}}\left(y\right)-{\mathit{P}}_{\mathbf{FOV}}\left(y\right)\right]}^2} \le R_1\\ R_1=H\bullet \tan \left(\frac{{\theta }_{\mathrm{FOV}}}{2}\right)+R_{\mathrm{tele}}+\frac{\left[{\mathit{P}}_{\mathit{i}}\left(z\right)-H\right]\bullet \sin \left(\frac{{\theta }_{\mathrm{FOV}}}{2}\right)}{\sqrt{n_w^2-{\sin }^2\left(\frac{{\theta }_{\mathrm{FOV}}}{2}\right)}}\end{array},\end{array}$$(14)

where P_FOV is the center of the FOV section, and the R_tele is the radius of the telescope. For vertical incidence, P_FOV(x) and P_FOV(y) are equal to x_st and y_st. For the oblique incidence, P_FOV(x) and P_FOV(y) are related to P_FOV(z). Figure 4a describes the vertical detection situation when the lidar system is at P_st = [0, 0, 0]. If the photon interacts with the target, P_i is replaced with P_touch to determine whether the interaction position satisfies the detection condition. The signal contribution of jth photon at position P_i is estimated by [46]$$\begin{array}{c}N\left(i,j\right)=∆\mathrm{\Omega }\tilde{\beta }\left({\theta }_p\right)T_s{W}_i\exp \left(-{\tau }_i\right),\end{array}$$(15)

Fig. 4 Geometry of photon detection of MCT model. (a) Estimation of water signal (green dotted line) and target reflection signal (red dotted line). (b) The shaded area under the target and non-detectable water signal (gray dotted line).

where $\tilde{\beta }\left({\theta }_p\right)$ represents the probability that a photon will be scattered from the current direction to the receiving direction. If the photon interacts with the target, $\tilde{\beta }\left({\theta }_p\right)$ is the Lambertian phase function. Otherwise, $\tilde{\beta }\left({\theta }_p\right)$ is the phase function of mixed water in equation (8). T_s is the sea surface transmission, W_i is the present weight of the photon. The solid angle of the detector ∆Ω can be calculated by the detector’s image position P_img = [x_st, y_st, −n_wH + H]. τ_i is the optical depth from position P_i to lidar position P_st.

However, a special scenario that needs to be considered is when the photon moves underneath the target. As shown in Figure 4, for the vertical detection, the covered extent can be described approximately as$$\begin{array}{c}\{\begin{array}{c}x\in \left[\frac{-L_C}{2},\frac{L_C}{2}\right]\\ y\in \left[{-R}_2,R_2\right]\\ z>H+d_w\end{array},\end{array}$$(16)

where R₂ = d_imgtan(θ_cov/2), d_img = z + H(n_w − 1) is distance from the detector’s image position P_img to the circle center P_C, θ_cov=2sin⁻¹(R_C·d_img). When the photon enters this area, the detectable trajectory is blocked by the target, and the photon can no longer make contributions to the signal.

2.5 Interaction with the target

Noise is an influential factor of detection, which mainly includes background signal of solar radiance and detector noise. N_b represents the background signal and is given by [70]$$\begin{array}{c}{N}_b=I_{b}A∆\lambda \Delta t{\eta }_{\mathrm{oe}}{\eta }_{\mathrm{qe}}\frac{\lambda \pi {\theta }_{\mathrm{FOV}}^2}{4h{c}_v},\end{array}$$(17)

where I_b is the solar spectral radiance reflected from atmosphere and the water’s surface, A is receiver aperture area, ∆λ is the bandwidth of the detector’s filter, ∆t is the sample time, η_qe is the quantum efficiency of the detector, η_oe is the optical efficiency of the lidar system, λ is the laser wavelength, c_v is the speed of light in vacuum, h is the Planck constant. The signal-to-noise ratio of the target signal can be calculated by [70]$$\begin{array}{c}\{\begin{array}{c}\delta N=\sqrt{N_{\mathrm{st}}+N_{\mathrm{sw}}+N_b+N_d}\\ \mathrm{SNR}=\frac{\sqrt{m}\left(N_{\mathrm{st}}-N_{\mathrm{sw}}\right)}{\delta N}\end{array},\end{array}$$(18)

where N_st is the lidar returned photoelectrons with target detection, N_sw is the lidar returned photoelectrons without target (only water detection), N_d is the dark count signal of the detector, and m is the number of the laser shots that need to be integrated.

3 Scanning oceanic lidar system

The airborne linear scanning oceanic lidar system (SOL) for underwater target detection is shown in Figure 5. We chose a 532.1 nm Nd:YAG laser as the light source for strong water penetration. The laser operates reliably over a temperature range of 0–40 °C. Under laboratory conditions at 25 °C, both the output power and pulse energy exhibit stability better than 2%, effectively reducing variations in pulse width and beam quality induced by environmental temperature fluctuations. The photomultiplier tubes (PMT) from the H10720-01 series manufactured by Hamamatsu Photonics are chosen as the detector for good sensitivity and wide dynamic range. The PMT exhibits high cathode radiant sensitivity in the green spectral band and operates reliably over a temperature range of 5–50 °C. Owing to its low dark current, the detector is capable of responding to photon-level signals. The hardware parameters are shown in Table. 1.

Fig. 5 (a) The schematic diagram of SOL system, the laser beam path is in green. Ap is short for aperture stop, CL is short for collimating lens, BF is short for band filter, FL is short for focusing lens, and PBS is short for polarizing beam splitter. (b) SOL is deployed on the UAV.

In the receiving subsystem, a Fresnel lens is selected as the telescope. The telescope has a rear focal length of 60.5 mm and an optical transmittance of up to 96%. The telescope collects the echoes and focuses them onto the rear focal plane. After aperture -stop limitation, the echoes pass through a collimating lens and a narrowband filter to enhance beam collimation and spectral selectivity, suppressing interference from other wavelengths. The conditioned echoes are then separated by a polarizing beam splitter, refocused by focusing lenses, and detected by PMTs, enabling high-sensitivity detection of the return signals.

In the scanning subsystem, a paraxial structure is chosen with the emission subsystem and the receiver subsystem separated at both ends to implement mechanical synchronous scanning. Besides, a reflective photoelectric incremental encoder is equipped to obtain angle information during scanning. It continuously outputs high-resolution TTL pulse signals in A, B, and Z phases. One full rotation corresponds to 163480 pulse signals, yielding a minimum angular resolution of 7.9 arcseconds, which meets the angular accuracy requirements of the UAV-borne lidar system.

4 Model validation

We first simulated water signals and compared them with lidar equation results for model theoretical validation, then simulated the underwater target signals and compared them with the measured profiles to verify the function of the photon’s interaction with target in the model.

4.1 Comparison with lidar equation

The oceanic lidar equation has been broadly applied to validate various lidar signal simulators [47, 49]. The depth-dependent lidar return signal of the lidar equation is [42, 43]$$\begin{array}{c}{N}_s\left(z_w\right)=\frac{E}{\mathrm{h\nu }}\bullet \frac{A}{{\left(n_{w}H+z_w\right)}^2}\bullet O{T}_a^2{T}_s^2{c}_v\Delta t_p{\eta }_{\mathrm{oe}}{\times \eta }_{\mathrm{qe}}{\beta }_{\pi }\left(z_w\right)\exp \left(-2{\int }_0^{z_w}{K}_{\mathrm{lidar}}\left(z^{\text{'}}\right)d{z}^{\text{'}}\right),\end{array}$$(19)

where z_w =z − H is water depth in this study, Ns is the received signal, E is the laser energy for one pulse, O is the overlap factor, ∆t_p is the pulse width, T_a is the transmission of atmosphere, v is the frequency of the laser, β_π(z) is the volume scattering coefficient calculated by (6) at the scattering angle of π. K_lidar is the lidar attenuation coefficient which is calculated by [43]$$\begin{array}{c}{K}_{\mathrm{lidar}}\left(z_w\right)=K_d\left(z_w\right)+\left[c_{\mathrm{total}}\left(z_w\right)-K_d\left(z_w\right)\right]\\ \exp \left[-0.85c_{\mathrm{total}}\left(z_w\right)D\left(z_w\right)\right],\end{array}$$(20)

where c_total(z_w) is the extinction coefficient, D_w(z_w) is the lidar spot diameter at different water depth, K_d(z_w) is the diffuse attenuation coefficient.

Because of the limitation of lidar equation in dealing with multiple scattering, the water environment has to be as clean as possible to reduce multiple scattering. Referring to Morel’s research, we choose the simulation seawater from the South Pacific anticyclonic gyre (SPG) which is oligotrophic and the “clearest” [57]. Figure 6a shows the chlorophyll concentration profile of SPG water, and the calculated inherent optical properties based on the chlorophyll concentration [43, 57–62].

Fig. 6 (a) Chlorophyll concentration profile of SPG water, (b) comparison between MC simulation signal and lidar equation signal.

The relative error used here, which describes how closely the simulated signals match the lidar equation signal, is as follows$$\begin{array}{c}\delta =\frac{\left|X_s-X_R\right|}{X_R}\times 100\%,\end{array}$$(21)

where X_S is the simulated value, and X_R is the result of lidar equation. We take the mean relative error at 100 m intervals for statistics. We set the number of photons to 10⁷ and take the average results of five simulations for comparison. For simplicity, the MC simulation result and the lidar equation result are normalized based on their respective sea surface signals.

Figure 6b compares the normalized lidar return signals from the lidar equation and the MCT model. A strong agreement between the two results is shown over most depth ranges, indicating the MCT model can accurately reproduce the expected lidar signal attenuation behavior. Slight deviations appear at the deeper water layers, which can be attributed to the cumulative effects of multiple scattering. These effects are inherently captured by the MCT model but simplified in the analytical lidar equation. These deviations can be further characterized through quantitative analysis. The mean relative error is less than 5% within a depth of 50 m and increases to approximately 19.4% at a depth of 100 m. Overall, these results demonstrate that the proposed MCT model achieves high accuracy while providing a more comprehensive and physically realistic description of underwater lidar signal propagation.

4.2 Experiment validation

The measurement experiments of SOL were conducted in the Thousand-Island Lake, Zhejiang province, China. The mean depth of the Thousand-Island Lake is about 34 m, and the maximum depth is more than 100 m. We chose the area near the Gangkou Bridge between two groups of mountains as the experimental area (Fig. 7a). On this occasion, the wind disturbance can be drastically avoided, and the water surface is wide enough. The inherent optical properties of water were measured by Spectral Absorption and Attenuation Sensor ac-s (WET Labs Inc.), and the results are shown in Figure 7b. The green metal barrel was selected as the underwater target. Its length is 0.8 m and its radius is 0.3 m. Its surface has a strong reflection and the reflectivity is about 0.5. As shown in Figure 7c, the ball float and the clump weight were used to fix the metal barrel at a water depth of about 7 m. The lidar SOL was mounted on the UAV DJI Flycart 30. The FPV camera on the UAV can help monitor the surface and assist lidar detection. The direction of flight was perpendicular to the long side of the target (Fig. 7d). Case 1 and case 2 each represent a group of laser beams making up a single linear scanning.

Fig. 7 Field experiment of SOL, (a) experimental area in Thousand-Island Lake, (b) the inherent optical properties of experimental water, (c) underwater placement of the target, (d) geometry of the SOL flying over the underwater target.

The SOL measurements were processed through a series of corrections. In the data pre-processing phase, a low-complexity method based on linear approximation of the leading edge (LLE) was used to extract the water surface precisely, then the slope distances both in the air and water could be calculated as well [71, 72]. Next, according to the angle information of the rotating motor, the slope distances were converted to vertical distances. The one-dimensional signals could be put together as two-dimensional images. Finally, the intensity of the lidar’s output was converted to the number of photoelectrons by calibrating the single-photon response of the detector in the lab.

Figures 8a and 8b show two typical scanning images, and the relevant observed geometries of the UAV that flew over the target are shown in Figure 7d, marked as 1 and 2. Figure 8a was obtained near the top of the target center (case 1), while Figure 8b was obtained from the top edge of the target (case 2). Both of the results had strong signals at a water depth of about 7 m, which revealed the good detection capabilities of SOL. The intensity of the target signals in case 2 was declining and became difficult to distinguish. We also found the signals from the shallow water were so strong that the value was beyond the maximum the data acquisition unit could record. It was the reason that signals from the water surface to about 2 m depth were almost the same.

Fig. 8 Comparing simulated results with scanned profiles of SOL, (a) and (b) are SOL’s measurements labeled laser beam of 1 and 2 in Figures 6d, 6c and 6d are corresponding simulation results. The color bar represents the number of simulated photoelectrons.

It is apparent to see the target’s signals at a depth of 7.0 m in the simulated results of Figures 8c and 8d, which is very consistent with the measured results in Figures 8a and 8b. The simulated target’s signals scanned in case 1 are stronger than in case 2. Other than the consistency of the target’s position, the quantity of simulated photoelectrons was consistent with the order of magnitude of measured photons. All of the photoelectron numbers first decreased from 10³ at the water surface to about 10^2.5 at a depth of 3.0 m, then decreased to 10¹ at 7.0 m. The numbers of photoelectrons from the target were around 10² in case 1 and around 10^1.5 in case 2. To sum up, the above comparisons demonstrated the correctness of the MCT model.

5 Analysis of the detection capability of SOL

Using the MCT model, the detection capabilities of SOL are further explored. The simulation experiments are carried out under two typical kinds of seawater: Jerlov II and Jerlov 3C [73, 74]. IOPs shown in Figures 9a and 9b are from WOOD (World-wide Ocean Optics Database). The lidar system parameters and the size of the target are the same as in Section 4.2.

Fig. 9 IOPs of Jerlov II (a), Jerlov 3C (b).

To evaluate the detection differences between two kinds of water, we simulate the signals generated by the target placed at different z-coordinates in case 1, as shown in Figures 10 and 11. The scanning horizontal resolution is 0.4 m. The simulated results show that, when the depth of the target is smaller than 25.0 m in Jerlov II water and 9.0 m in Jerlov 3C water, it is easy to recognize the target. When the z-coordinate of the target is 30.0 m in Jerlov II water and 11.0 m in Jerlov 3C water, the target is barely possible to recognize. Besides, the thickness of the target’s signal, which means its depth range, is greater in Jerlov 3C. This might be relevant to the multiple scattering.

Fig. 10 Simulated results of target detection at different depths in Jerlov II water. The color bar represents the number of simulated photoelectrons.

Fig. 11 Simulated results of target detection at different depths in Jerlov 3C water. The color bar represents number of simulated photoelectrons.

Generally speaking, the target recognition of a two-dimensional image is more difficult than that of a single profile. This is because in a two-dimensional image, there must be continuous target features (at least two continuous profiles) to be convincing. As shown in Figures 10 and 11, the strongest target signals exist near the center of it (x = 0), and decrease with horizontal distance from the center. Based on the results of Figures 10 and 11, we simulate the detection capability of a single profile of multiple z-coordinates of the target in case 1. We use N_st−N_sw to acquire the signals that are only reflected by the target. The feature peaks of N_st−N_sw signals with different depths are extracted and plotted as Figure 12a. The signal-to-noise ratios of these feature peaks are calculated and shown as Figure 12b. It can be seen that the photoelectron numbers in both water environments go down exponentially with depth, and in Jerlov 3C water, they decrease much faster because of the bigger attenuation coefficients. The same trend appears in the distribution of target signal SNR with respect to depth, as shown in Figure 12b.

Fig. 12 Analysis of simulation results for target detection at different depths in Jerlov II water and Jerlov 3C water, (a) simulated numbers of photoelectrons, (b) SNR of target signals, (c) extended detection range of different depths, (d) horizontal resolution of the scanning in x-axis of different depths.

It is widely known that SNR can be improved by integrating several profiles. To improve the SNR of N_st−N_sw, the horizontal sampling interval should be reduced so that more laser shots are able to encounter the target. We have already discovered that the horizontal detectable length range of target is beyond its actual length [xϵ (−0.4 m, 0.4 m)] (Figs. 10 and 11). Therefore, after counting the number of photoelectrons with horizontal distribution at the same depth of target, the definition of extended detection range is made as the x-range when the photoelectron numbers with target reduce to half of that at the center top (x = 0). The statistical results of the extended detection range are displayed in Figure 12c. Assume the SNR does not change within this range, and set SNR = 3 as the standard. Then the relationships between max detectable depth and horizontal scanning resolution of Jerlov II and Jerlov 3C water are shown in Figure 12d by interpolating at the depth resolution of 0.1 m. Only one single profile is needed when the target depth is shallower than 25.0 m in Jerlov II water and 9.2 m in Jerlov 3C water, so the max point resolutions are equal to the extended detection ranges. It should be noted that the phenomenon of photoelectron numbers of the target dropping to half does not mean the target cannot be recognized. Recognition still depends on its true value and SNR (Fig. 10a). As the depths of target add up, more and more profiles have to be integrated, so the point resolutions are equal to the corresponding extended detection range divided by the number of integrated profiles. If the SNR reduces to a very small value, like lower than 1, it will be better to level up the number of integrated four times as calculated by equation (18).

The preceding analysis reveals that, to achieve a certain target SNR, deeper targets require the integration of more lidar profiles to enhance signal quality. This requires an increase in horizontal sampling density to ensure that more laser beams effectively interact with the target. While this approach improves the detectability of deep targets, it simultaneously reduces the area that can be covered within a given time, thereby decreasing the overall detection efficiency. This compromise carries important scientific implications for the design and operation of oceanic scanning lidar systems for target detection. Firstly, it quantitatively links detection depth to the required lidar signal profiles integration, and offers a theoretical basis for estimating detection costs at varying depths. Secondly, it highlights the inherent constraint between detection depth and operational efficiency and emphasizes the balance of spatial resolution, scanning speed, and signal quality in system configuration. This understanding provides valuable guidance for making acquisition strategies to different targets and water environments.

6 Conclusion

Oceanic lidar is significant for detecting deeper water. In this study, we proposed a semi-analytic MC model of underwater target detection (MCT) and developed an airborne linear scanning oceanic lidar system (SOL) for model validation. The comparison is made between lidar equation signals and MC simulated signals to validate the correctness of our model. It turns out the mean relative error is less than 5% within a depth of 50 m and is about 19.4% within a depth of 100 m. By comparing the simulated two-dimensional scanning signals with SOL’s field target measurements, the measured results show the consistency of the target’s position and the order of magnitude of photoelectrons. Then a series of simulations in Jerlov II water and Jerlov 3C water are made to study the detection capabilities of SOL. It turns out that, set SNR = 3 as the detection standard, the maximum target detection depths are 25.0 m in Jerlov II water and 9.2 m for a single profile. The definition of extended detection range is made, and how the maximum horizontal scanning resolution varies with depth is studied as well. Their relationships are conducive to balancing detection capabilities and scanning efficiency, providing a good reference for practical detection scenes. By exploring a range of representative simulation conditions, including different water optical types, target depths, and scanning parameters, the proposed MCT model is shown to be robust and effective. In addition to the cases presented in this work, the model framework is readily extensible to other scenarios involving variations in lidar system parameters and target-related properties, such as target size, shape and reflectivity.

There are also limitations to the MCT, such as the floating ball interference. In practical field measurements, floating balls can generate strong reflective signals at the sea surface, introducing additional non-target returns and locally increasing background levels. Such surface-related reflections may reduce the contrast of target echoes under certain conditions. In addition, other detector-related noise such as the after pulse, manifests primarily as weak delayed components in the received waveforms, which may cause a slight overestimation of deep-water signals. These effects will be further investigated in future studies. Furthermore, a more comprehensive parametric study about diverse water optical properties and lidar system parameters will be essential to construct generalized models, which can improve the transferability and applicability of the simulated findings in the broader practical scenarios.

Funding

This work is supported by the National Key Research and Development Program of China (Grant No. 2022YFB3901701).

Conflicts of interest

The authors declare that there are no financial interests, commercial affiliations, or other potential conflicts of interest that could have influenced the objectivity of this research or the writing of this paper.

Data availability statement

The data underlying the results presented in this paper is not publicly available at this time but may be provided upon reasonable request to the corresponding author.

Author contribution statement

All authors have reviewed, discussed, and agreed to their individual contributions. The specific situation is as follows: Methodology, Xinke Hao and Yan He; Data Curation, Xinke Hao, Yingjie Ruan and Hui Qi; Writing – Original Draft Preparation, Xinke Hao and Huixin He; Writing – Review and Editing, Huixin He and Yan He; Supervision, Deliang Lv and Guangxiu Xu; Project Administration, Yan He and Junwu Tang; Funding Acquisition, Yan He.

References

All Tables

All Figures

	Fig. 1 Flow chart of MCT simulation model.
In the text

	Fig. 2 (a) Global coordinate system of the simulation model, and (b) sampling the point from the Gaussian spot.
In the text

	Fig. 3 Geometry of the photon colliding with the cylinder target, the left sub-graph is projection of the interaction process.
In the text

	Fig. 4 Geometry of photon detection of MCT model. (a) Estimation of water signal (green dotted line) and target reflection signal (red dotted line). (b) The shaded area under the target and non-detectable water signal (gray dotted line).
In the text

	Fig. 5 (a) The schematic diagram of SOL system, the laser beam path is in green. Ap is short for aperture stop, CL is short for collimating lens, BF is short for band filter, FL is short for focusing lens, and PBS is short for polarizing beam splitter. (b) SOL is deployed on the UAV.
In the text

	Fig. 6 (a) Chlorophyll concentration profile of SPG water, (b) comparison between MC simulation signal and lidar equation signal.
In the text

	Fig. 7 Field experiment of SOL, (a) experimental area in Thousand-Island Lake, (b) the inherent optical properties of experimental water, (c) underwater placement of the target, (d) geometry of the SOL flying over the underwater target.
In the text

	Fig. 8 Comparing simulated results with scanned profiles of SOL, (a) and (b) are SOL’s measurements labeled laser beam of 1 and 2 in Figures 6d, 6c and 6d are corresponding simulation results. The color bar represents the number of simulated photoelectrons.
In the text

	Fig. 9 IOPs of Jerlov II (a), Jerlov 3C (b).
In the text

	Fig. 10 Simulated results of target detection at different depths in Jerlov II water. The color bar represents the number of simulated photoelectrons.
In the text

	Fig. 11 Simulated results of target detection at different depths in Jerlov 3C water. The color bar represents number of simulated photoelectrons.
In the text

	Fig. 12 Analysis of simulation results for target detection at different depths in Jerlov II water and Jerlov 3C water, (a) simulated numbers of photoelectrons, (b) SNR of target signals, (c) extended detection range of different depths, (d) horizontal resolution of the scanning in x-axis of different depths.
In the text