A Novel De-Noising Method for Improving the Performance of Full-Waveform LiDAR Using Differential Optical Path

A novel de-noising method for improving the performance of full-waveform light detection and ranging (LiDAR) based on differential optical path is proposed, and the mathematical models of this method are developed and verified. Backscattered full-waveform signal (BFWS) is detected by two avalanche photodiodes placed before and after the focus of the focusing lens. On the basis of the proposed method, some simulations are carried out and conclusions are achieved. (1) Background noise can be suppressed effectively and peak points of the BFWS are transformed into negative-going zero-crossing points as stop timing moments. (2) The relative increment percentage of the signal-to-noise ratio based on the proposed method first dramatically increases with the increase of the distance, and then the improvement gets smaller by increasing the distance. (3) The differential Gaussian fitting with the Levenberg-Marquardt algorithm is applied, and the results show that it can decompose the BFWS with high accuracy. (4) The differential distance should not be larger than c/2 × τrmin, and two variable gain amplifiers can eliminate the inconsistency of two differential beams. The results are beneficial for designing a better performance full-waveform LiDAR.


Introduction
Light Detection and Ranging (LiDAR) is an active, remote-sensing system that provides direct range measurement and is capable of collecting three-dimensional (3D) spatial information [1,2]. The system emits short laser pulses with high frequency to illuminate the object surface and then employs photodiode detectors to record backscattered waveform signals [3][4][5]. LiDAR has received much attention in recent years due to its simplicity, high accuracy, and utility in many applications, such as 3D city modeling, cartography, forest inventories, target recognition, and the digital terrain model [6][7][8][9]. To the best of our knowledge, LiDAR can be divided into two categories, i.e., discrete-echo LiDAR and full-waveform LiDAR [10,11]. A discrete-echo LiDAR records only a few discrete-echo backscattered waveform signals for each transmitted laser pulse and provides only 3D coordinates and range information about objects [8,12]. Unlike the discrete-echo LiDAR, a full-waveform LiDAR can record the entire backscattered full-waveform signal (BFWS) for each transmitted laser pulse

Principle
The principle of the full-waveform LiDAR based on the differential optical path is shown in Figure 1.
First, a field programmable gate array (FPGA) generates a trigger signal for a laser to emit a laser pulse, and the laser pulse is collimated by a transmitting lens.
Second, the transmitted laser pulse is divided into two beams by beam splitter 1. (1) One beam is focused by a convergent lens and is detected directly by a photo detector [30]. Because the beam undergoes no noise interference and it keeps the original stand waveform, it is suitable to regard the peak position of the laser pulse as the start timing moment for a timer in FPGA. (2) The other beam is projected into a scene for illuminating objects. Three objects are assumed to exist in the scenario, and the distances between the three objects and the laser are R 1 , R 2 , and R 3, respectively. The BFWS-the sum of the backscattered sub-waveform signal (BSWS)-is generated after the interactions of the transmitted laser pulse with each encountered object.
Remote Sens. 2017, 9,1109 3 of 16 Third, the BFWS is reflected by the beam splitter 1 and focused by a focusing lens. Two APDs are placed before and after the focus of the focusing lens (offset distance is L). The beam splitter 2 divides the BFWS into two beams and the two beams impinge on the two APDs, respectively. The two electrical signals (P r1 and P r2 ) of the BFWS are sent to a subtraction circuit (SC). After the subtraction operation of the SC, a differential BFWS (P rd ) is obtained and sampled by an analog-to-digital converter (ADC). The differential BFWS has some negative-going zero-crossing points (NGZCPs) and it is set as the stop timing moments for the timer in the FPGA.
Finally, the other parameters of the different objects, including the amplitude, the position traveling time, and the standard deviation, are obtained by analyzing the differential BFWS. The time of flights are determined between the start and stop timing moments. The range information about the objects can be obtained using the time of flights. Unlike our previous study [21], besides the range information, other information related to geometry and radiative properties of the illuminated objects can be achieved by the other parameters through decomposing the differential BFWS.
position is set as the start timing moment. In terms of the stop timing moments, the traditional method uses only a single APD for detecting stop signals. The stop timing moments of different objects are the peak points (PPs) of the BFWS, as shown in Figure 2a. However, the peak positions are difficult to discriminate for the timer, because the temporal change rates near the peak positions are small, which increases the difficulty in obtaining the peak position. Unlike the traditional method, the differential optical path method employs two APDs to receive the BFWS. These two APDs are placed before and after the focus of the focusing lens, respectively. Compared with the peak discriminator of the stop timing moments based on the traditional method, the stop timing moments for the timer are changed into the NGZCPs of the differential BFWS based on differential optical path method, which is shown in Figure 2b. We can see that the temporal change rates of the differential BFWS at the NGZCPs in Figure 2b are obviously higher than those at the PPs of the BFWS in Figure  2a. Therefore, the stop timing moments can be easily detected by the timer through the differential optical path method. From the principle of the proposed method, we can see that the amplitudes of the two BFWSs impinging on by the two APDs (Pr1 and Pr2) are reduced by half compared with those of a single APD. Hence, in order to ensure the proposed method works properly, the amplitudes of the two BFWSs should be larger than the minimum input power of the APD. Although the amplitudes of the two BFWSs are reduced by half, the background noise can be suppressed and the SNR can be enhanced by subtracting the two BFWSs. Therefore, the range information and other information about the objects based on the proposed method are more accurate than that of the traditional method. Figure 1. The principle of the full-waveform light detection and ranging (LiDAR) system based on the differential optical path. Figure 1. The principle of the full-waveform light detection and ranging (LiDAR) system based on the differential optical path.
The range information is an important parameter obtained from the BFWS and is determined by the start and stop timing moments. The start timing moment comes from the photo detector, and the peak position is easy to detect because of the high intensity and low background power, so the peak position is set as the start timing moment. In terms of the stop timing moments, the traditional method uses only a single APD for detecting stop signals. The stop timing moments of different objects are the peak points (PPs) of the BFWS, as shown in Figure 2a. However, the peak positions are difficult to discriminate for the timer, because the temporal change rates near the peak positions are small, which increases the difficulty in obtaining the peak position. Unlike the traditional method, the differential optical path method employs two APDs to receive the BFWS. These two APDs are placed before and after the focus of the focusing lens, respectively. Compared with the peak discriminator of the stop timing moments based on the traditional method, the stop timing moments for the timer are changed into the NGZCPs of the differential BFWS based on differential optical path method, which is shown in Figure 2b. We can see that the temporal change rates of the differential BFWS at the NGZCPs in Figure 2b are obviously higher than those at the PPs of the BFWS in Figure 2a. Therefore, the stop timing moments can be easily detected by the timer through the differential optical path method. From the principle of the proposed method, we can see that the amplitudes of the two BFWSs impinging on by the two APDs (P r1 and P r2 ) are reduced by half compared with those of a single APD. Hence, in order to ensure the proposed method works properly, the amplitudes of the two BFWSs should be larger than the minimum input power of the APD. Although the amplitudes of the two BFWSs are Remote Sens. 2017, 9, 1109 4 of 16 reduced by half, the background noise can be suppressed and the SNR can be enhanced by subtracting the two BFWSs. Therefore, the range information and other information about the objects based on the proposed method are more accurate than that of the traditional method.

Analysis of the Differential BFWS
Different from our previous study [21], the fundamental expression to describe the differential BFWS is based on the LiDAR equation. The power of the transmitted laser pulse is supposed to be a temporal function of Gaussian model and it can be written as [31] ( ) where Et is the original pulse energy, and τ is the transmitting pulse width. The pulse laser transmits a narrow laser beam with a certain divergence angle toward the object, and the power impinging on the APD is written as [3] where R is the distance between the object and laser, βt is the transmitter beam divergence, ρ is the reflectivity, As is the receiving area of the object, Ω is the solid angle, and Dr is the aperture diameter of the receiver optics system. Equation (2) can also be rewritten as the following LiDAR equation [32] where σ is the backscatter cross-section, which represents the character of the object, such as the reflectivity and the directionality of scattering. From a practical point of view, additional power losses in the instrument and atmosphere must be considered. Therefore, the ultimate LiDAR equation, i.e., the BSWS, is [12,32] where ηsys is the system transmission factor, and ηatm is the atmospheric transmission factor. As mentioned above, the traditional method employs a single APD to receive the BFWS. Therefore, if N distinct objects exist within the travel path of the laser pulse in the scene, then the

Analysis of the Differential BFWS
Different from our previous study [21], the fundamental expression to describe the differential BFWS is based on the LiDAR equation. The power of the transmitted laser pulse is supposed to be a temporal function of Gaussian model and it can be written as [31] where E t is the original pulse energy, and τ is the transmitting pulse width. The pulse laser transmits a narrow laser beam with a certain divergence angle toward the object, and the power impinging on the APD is written as [3] where R is the distance between the object and laser, β t is the transmitter beam divergence, ρ is the reflectivity, A s is the receiving area of the object, Ω is the solid angle, and D r is the aperture diameter of the receiver optics system. Equation (2) can also be rewritten as the following LiDAR equation [32] where σ is the backscatter cross-section, which represents the character of the object, such as the reflectivity and the directionality of scattering. From a practical point of view, additional power losses in the instrument and atmosphere must be considered. Therefore, the ultimate LiDAR equation, i.e., the BSWS, is [12,32] Remote Sens. 2017, 9, 1109 5 of 16 where η sys is the system transmission factor, and η atm is the atmospheric transmission factor. As mentioned above, the traditional method employs a single APD to receive the BFWS. Therefore, if N distinct objects exist within the travel path of the laser pulse in the scene, then the expression of the BFWS, i.e., the sum of the BSWSs of each distinct object can be written as [3] where R i is the distance between the i-th object and the laser, τ ri is the received pulse width, c is the light speed, θ i is the tilt angle of the i-th object between normal vector and optical axis, W 0 is the waist radius of the laser, and W(R i ) is the beam radius at the R i distance.
According to the principle based on the differential optical-path as shown in Figure 1, the BFWSs of the APD 1 and the APD 2 are written as Equation (6) is expressed under ideal condition, i.e., the noise is ignored. The BGN is an external noise of the APD and is a major noise source especially in the bright sunlight measurement environment. Therefore, it should be taken into consideration. The power of the received background solar illumination impinging on the APD can be written as [20] P Bi = ρ i h sum T r A r sin(FOV/2) 2 ∆λ, where h sum is the background solar irradiance, T r is the transmission of the receiver, FOV is the field of view of the receiver optics system, A r is the area of the receiver, and ∆λ is the optical bandwidth. Equation (7) shows that the BGN is a constant, i.e., it does not vary with time. Moreover, according to the previous studies [33,34], in LiDAR measurement systems, the level of background noise generated by the atmosphere can be treated as constant. Therefore, the two BFWSs impinging on the APD 1 and APD 2 with the BGN are written as After the subtraction operation between these two BFWSs detected by the APD 1 and the APD 2, the differential BFWS is obtained and it can be written as Equation (9) shows that the BGN is suppressed effectively. According to our previous study [21], the time of flight of each object using the differential optical path method also equals to 2R i /c. In other words, the range information is unaffected by the proposed method. Moreover, compared with traditional peak discriminator, the PP is changed into the NGZCP, which is more easily detected by the Remote Sens. 2017, 9, 1109 6 of 16 timer than that of the PP. Therefore, the differential optical path method can more easily obtain range information and other information than the traditional method.

SNR Analysis
The SNR is a synthetical parameter to evaluate the quality of the full-waveform data [1]. The SNR expression of a single APD based on the traditional method is written as [21,33,35] where P sig is the detected signal power of the APD, P th is the mean-squared thermal-noise power, P a is the mean-squared noise power added by the electronic amplifier, P dark is the mean-squared dark-current-noise power of the APD, P shot is the mean-squared signal shot noise power, and P back is the output power of background noise by the APD. The expressions of the above terms based on the traditional method are described as follows where I is the detected photocurrent, ρ D = η D e/hf is the current responsivity of the APD, η D is the quantum efficiency of the APD, e is the electron charge, h is Planck's constant, f = c/λ is the frequency of light, R L is the effective load resistance of the APD, M is the current gain of the APD, ρ D is the responsivity of the APD, k is Boltzmann's constant, B is the electrical bandwidth of the system, i.e., B = 1/(2τ), T is the temperature in Kelvin, T a is the effective noise temperature, I dark is the dark current, and F ex is the excess-noise factor. By substituting Equation (11) in to Equation (10), we obtain the SNR expression of the traditional method as which can be rewritten as Employing the relationship between the responsivity and quantum efficiency of the APD, i.e., ρ D = η D e/hf, the SNR of the traditional method also be written as Compared with the traditional method, the proposed method employs two APDs to receive the BFWS. Because the background noise P back is constant, it can be suppressed by subtraction operation of the two BFWSs detected by two APDs. However, the other four noise are not constants; in fact, Remote Sens. 2017, 9, 1109 7 of 16 they are random numbers whose possibility distribution functions obey Poisson distributions [36]. Therefore, the SNR expression of the proposed method with two APDs can be achieved by where ρ D1 and ρ D2 are the responsivities of the two APDs, R L1 and R L2 are the effective load resistances of the two APDs, M 1 and M 2 are the current gains of the two APDs, T 1 and T 2 are the temperatures in Kelvin of the two APDs, T a1 and T a2 are the effective noise temperatures of the two APDs, F ex1 and F ex2 are the excess noise factors of the two APDs, and I dark1 and I dark2 are the dark currents of the two APDs.

Waveform Decomposition and Differential Gaussian Fitting
As mentioned above, compared with the discrete-echo LiDAR, the full-waveform LiDAR can provide more potential ability in extracting additional parameters and deriving properties of the objects from the BFWS. Decomposition is the core process for deriving the valuable parameters [37]. Gaussian fitting is one of the most commonly adopted methods for decomposing full-waveform data [38,39]. The Gaussian fitting is based on the assumption that the transmitted pulse is of a Gaussian type and the BFWS is composed of several single echoes that also are Gaussian types. Therefore, the Gaussian fitting expression of the BFWS based on the traditional method can be written as [12,32] where b is the noise level of the waveform, N is the peak number of the BFWS, i.e., the echo number, a i is the amplitude of the i-th echo, t i is the traveling time of the i-th echo, δ i is the half width of the i-th echo (standard deviation), and t is the traveling time.
The initial positions and echo number should be determined prior to the Gaussian fitting for the iteration process. In general, two conventional methods, i.e., the center of gravity and zero-crossing of the first derivative, can be used. Compared with conventional methods, the second derivative algorithm outperforms the two aforementioned traditional methods [40]. The second derivative of the BFWS is calculated as [37] where t i indicates an echo location of the BFWS, and ∆t is the time interval. In the second-derivative algorithm, a local minimum point is supposed to be the initial position of the echo and represents one echo number. Unlike the Gaussian fitting of the tradition method, the differential Gaussian fitting expression of the differential BFWS based on differential optical path method is written as Compared with Equations (5) and (18), we can see that the differential Gaussian fitting can preserve the characteristics of the BFWS, including the amplitude, the distance of the objects, and the pulse width. In terms of the amplitude, the amplitude of the differential Gaussian fitting is reduced to half of the BFWS due to the existing BS 2; therefore, it is needed to multiply the amplitude of the differential Gaussian fitting by two in order to correctly retrieve the parameters of the objects, including the backscatter cross-section. In terms of the distance of the objects, the NGZCP locations are the points where P rd (t) = 0 in Equation (18), which can be obtained by Then, we can achieve t = t i = 2R/c, which illustrates that the distance of the objects can be preserved by the differential Gaussian fitting. In terms of the pulse width, we can find that the standard deviation in Equation (18) is equal to the pulse width in Equation (5), and the tilt angle of the objects can be retrieved by the standard deviation.
The aforementioned method uses the second derivative of the BFWS to achieve echo number. On the contrary, the differential optical-path method regards the NGZCP number of differential BFWS as the echo number, which is written as If the sampling point P rd (t i ) of the differential BFWS satisfies Equation (19), then this sampling point is supposed to be the initial position of the echo and represents one echo number.
A nonlinear least-squares method with robust Levenberg-Marquardt (LM) technique is used to obtain the additional parameters (amplitude, position, and standard deviation) from the differential Gaussian fitting expression in Equation (18) [32]. The fitting quality is evaluated by a variable ξ, which is written as where C is a weight value and equals to 1/N, and ω is the desired accuracy determined by the end user.

Simulation Parameters and Model Verification
On the basis of the abovementioned analysis, the simulations of the BSWS of each object and the BFWS are carried out. Given that the BGN is one of the major noise sources, it is also simulated. Three objects with different reflectivities are selected and are positioned at different ranges to the laser. The set parameters of the pulse laser and the three objects are as shown in Table 1 [3, 16,20].
The BSWS and the BGN of each object based on the traditional method using the aforementioned parameters in Table 1 are shown in Figure 3. The Figure 3 shows that the power of BGN is higher than that of BSWS under the parameters, i.e., the BSWS is submerged under the BGN. Therefore, detecting the PPs of the BSWS is difficult. In order to solve this issue and test the model verification of the proposed method, the differential BSWS of each target is achieved using the differential optical-path method (the differential distance L = 0.03 m), and the corresponding results are shown in Figure 4a. From Figure 4a, we can see that the BGN is suppressed effectively and the PPs are transformed into NGZCPs. than that of BSWS under the parameters, i.e., the BSWS is submerged under the BGN. Therefore, detecting the PPs of the BSWS is difficult. In order to solve this issue and test the model verification of the proposed method, the differential BSWS of each target is achieved using the differential opticalpath method (the differential distance L = 0.03 m), and the corresponding results are shown in Figure  4a. From Figure 4a, we can see that the BGN is suppressed effectively and the PPs are transformed into NGZCPs. The differential BFWS, i.e., the sum of the differential BSWS of each object in Figure 4a is shown in Figure 4b. From Figure 4b, the number of the objects and the positions of the NGZCPs according to Equation (19) can be obtained. We can see that there are three NGZCPs that exist in the differential BFWS in Figure 4b. Therefore, the number of the objects is three, which equals to the object number in Table 1. According to the simulation parameters (the distances to the laser) in Table 1, the times of flights of the three objects are 3.3333, 3.3340, and 3.3353 μs, which are equal to the NGZCP positions of the differential BFWS, shown in Figure 4b. The aforementioned analysis shows that the BG is suppressed effectively and the time of flight, i.e., the range information of each object is unaffected by the proposed method.  Figure 4b. Therefore, the number of the objects is three, which equals to the object number in Table 1. According to the simulation parameters (the distances to the laser) in Table 1, the times of flights of the three objects are 3.3333, 3.3340, and 3.3353 μs, which are equal to the NGZCP positions of the differential BFWS, shown in Figure 4b. The aforementioned analysis shows that the BG is suppressed effectively and the time of flight, i.e., the range information of each object is unaffected by the proposed method.  The differential BFWS, i.e., the sum of the differential BSWS of each object in Figure 4a is shown in Figure 4b. From Figure 4b, the number of the objects and the positions of the NGZCPs according to Equation (19) can be obtained. We can see that there are three NGZCPs that exist in the differential BFWS in Figure 4b. Therefore, the number of the objects is three, which equals to the object number in Table 1. According to the simulation parameters (the distances to the laser) in Table 1, the times of flights of the three objects are 3.3333, 3.3340, and 3.3353 µs, which are equal to the NGZCP positions of the differential BFWS, shown in Figure 4b. The aforementioned analysis shows that the BG is suppressed effectively and the time of flight, i.e., the range information of each object is unaffected by the proposed method.

SNR Improvement
SNR is used to evaluate the quality the BFWS because it affects the measurement accuracy of the full-waveform LiDAR. Equations (14) and (15) indicate that the SNR is affected by many factors, such as dark current, temperature, current gain, and load resistance of the APD. Moreover, according to Equation (9), the power of the BFWS changes with the distance to laser and the power of the BGN changes the reflectivity of the objects. Therefore, considering so many factors at the same time is a challenge. Actually, the parameters of the two APDs and the pulse lasers can usually be determined when the full-waveform LiDAR system is already given. The SNR is related only to distance to the laser and reflectivity of the objects. To obtain the SNR, the reflectivity and the distance difference of the objects are fixed, shown in Table 1. The distance between the pulse laser and the first object varies from 300 m to 3000 m. Some typical parameters of the two APDs are set as follows: η D = η D1 = 0.9, η D2 = 0.8, e = 1.602 × 10 −19 C, h = 6.63 × 10 −34 J·s, and k = 1.38 × 10 −23 J/K [33].
The results of the SNR of the traditional method and the proposed method are shown in Figure 5a. The Figure 5a shows that: (1) The SNRs of the two methods decrease dramatically with the increase of the distance when the distance is shorter than 1000 m. However, the SNRs decrease slightly when the distance is shorter than 1000 m. (2) Based on the traditional method, the SNR decreases from 63 dB to 3 dB when the distance increases from 300 m to 1000 m and decreases from 3 dB to 0.009 dB when the distance increases from 1000 m to 3000 m. Based on the proposed method, the SNR decreases from 104 dB to 7 dB when the distance increases from 300 m to 1000 m, and decreases from 7 dB to 0.023 dB when the distance increases from 1000 m to 3000 m. (3) Compared with the traditional method, the proposed method can improve the SNR effectively at the same distance. For example, at a distance of 300 m, the SNR of the traditional method is 63 dB, whereas that of the proposed method is 104 dB. To clearly illustrate the SNR improvement by the proposed method, relative increment percentages of the SNR, i.e., ∆SNR = [(SNR p − SNR t )/SNR t ] × 100%, are calculated, where SNR p and SNR t are the SNRs of the proposed method and the traditional method, respectively. The results of the ∆SNR are shown in Figure 5b. The Figure 5b shows that the ∆SNR increases with the increase of the distance. When the distance increases from 300 m to 1000 m, the ∆SNR increases dramatically from 65 to 145. However, it increases slightly from 145 to 161 when the distance increases from 1000 m to 3000. The results show that the proposed method can improve the SNR, but the improvement gets smaller by increasing distance.
3 dB to 0.009 dB when the distance increases from 1000 m to 3000 m. Based on the proposed method, the SNR decreases from 104 dB to 7 dB when the distance increases from 300 m to 1000 m, and decreases from 7 dB to 0.023 dB when the distance increases from 1000 m to 3000 m. (3) Compared with the traditional method, the proposed method can improve the SNR effectively at the same distance. For example, at a distance of 300 m, the SNR of the traditional method is 63 dB, whereas that of the proposed method is 104 dB. To clearly illustrate the SNR improvement by the proposed method, relative increment percentages of the SNR, i.e., ΔSNR = [(SNRp − SNRt)/SNRt] × 100%, are calculated, where SNRp and SNRt are the SNRs of the proposed method and the traditional method, respectively. The results of the ΔSNR are shown in Figure 5b. The Figure 5b shows that the ΔSNR increases with the increase of the distance. When the distance increases from 300 m to 1000 m, the ΔSNR increases dramatically from 65 to 145. However, it increases slightly from 145 to 161 when the distance increases from 1000 m to 3000. The results show that the proposed method can improve the SNR, but the improvement gets smaller by increasing distance.

Waveform Decomposition and Differential Gaussian Fitting Accuracy
Waveform decomposition and Gaussian fitting are widely used methods in the full-waveform LiDAR and are the important processes for deriving the valuable parameters and features of the objects, such as amplitude, position (traveling time), and standard deviation from the differential BFWS. The differential BFWS should be decomposed to obtain these parameters. A differential Gaussian fitting with the LM technique for the proposed method is used, and the desired accuracy ω is set as 1 × 10 −16 . The real value of the differential BFWS is shown in Figure 4b. The differential

Waveform Decomposition and Differential Gaussian Fitting Accuracy
Waveform decomposition and Gaussian fitting are widely used methods in the full-waveform LiDAR and are the important processes for deriving the valuable parameters and features of the objects, such as amplitude, position (traveling time), and standard deviation from the differential BFWS. The differential BFWS should be decomposed to obtain these parameters. A differential Gaussian fitting with the LM technique for the proposed method is used, and the desired accuracy ω is set as 1 × 10 −16 . The real value of the differential BFWS is shown in Figure 4b. The differential Gaussian fitting results of the differential BFWS about the three objects and fitting accuracy are shown in Table 2. The results show that the relative errors of the amplitude of each object are 0.41%, 0.78%, and 0.29%, respectively; the relative errors of the position of each object are all 0%, which shows the proposed method can precisely detect the position of the objects; the relative errors of the standard deviation of each object are 0.07%, 0.10%, 0.01%, respectively; and the relative errors of the backscatter cross-section of each object are 0.51%, 0.89%, and 0.34%, respectively.
The differential Gaussian fitting curves of the BSWS of each object are shown in Figure 6a. The differential Gaussian fitting values (shown in Figure 6a) and the real values BFWS of each object (shown in Figure 4a) are compared. The absolute errors of these values are shown in Figure 6b. The maximum absolute errors of the three objects are 3.56 × 10 −12 , 2.80 × 10 −11 , and 2.08 × 10 −12 W, respectively. The minimum absolute errors of the three objects are −6.23 × 10 −12 , −1.99 × 10 −11 , and −4.58 × 10 −12 W, respectively. These results indicate that the differential Gaussian fitting with the LM technique is capable of decomposing the differential BFWS with high accuracy and can well preserve the characteristics of the objects.  The differential Gaussian fitting curves of the BSWS of each object are shown in Figure 6a. The differential Gaussian fitting values (shown in Figure 6a) and the real values BFWS of each object (shown in Figure 4a) are compared. The absolute errors of these values are shown in Figure 6b. The maximum absolute errors of the three objects are 3.56 × 10 −12 , 2.80 × 10 −11 , and 2.08 × 10 −12 W, respectively. The minimum absolute errors of the three objects are −6.23 × 10 −12 , −1.99 × 10 −11 , and −4.58 × 10 −12 W, respectively. These results indicate that the differential Gaussian fitting with the LM technique is capable of decomposing the differential BFWS with high accuracy and can well preserve the characteristics of the objects.

Discussion
Equation (9) shows that the expression of the differential BFWS based on the proposed method is affected by the differential distance (L) of the APD 1 and the APD 2. In other words, different differential distances exhibit different expressions of the differential BFWS. According to the

Discussion
Equation (9) shows that the expression of the differential BFWS based on the proposed method is affected by the differential distance (L) of the APD 1 and the APD 2. In other words, different differential distances exhibit different expressions of the differential BFWS. According to the principle of the Figure 1, an overlapping area should exist between the two BFWSs detected by the APD 1 and the APD 2 to ensure that the proposed system works properly. Meanwhile, given that the complexity of the proposed system, inconsistent splitting ratio of the BS3, and the individual difference of the two APDs (APD1 and APD2), it is challenging to maintain the consistent characteristics of the two beams split by the BS 2. Therefore, the inconsistency of the two beams should be eliminated prior to utilization.

Differential Distance Selection
An overlapping area should exist between the two BFWSs detected by the two APDs. Therefore, the differential time (2L/c) should not exceed the minimum value of the received pulse width (τ rmin ). Thus, the differential distance should not exceed c/2 × τ rmin , shown in Figure 7a. The differential BFWS of the two APDs is shown in Figure 7b when the differential distance is larger than c/2 × τ rmin . The Figure 7b shows that the number of the NGZCPs is two, i.e., the number of the object is two. In fact, three objects are included in the simulated scene. Therefore, a differential distance that is larger than c/2 × τ rmin may cause misjudgment of the number of the objects. Moreover, the differential BFWS is distorted, which increases the difficulty of wave decomposition and decreases the differential Gaussian fitting accuracy. In terms of the minimum value of the received pulse width τ rmin , it can be evaluated before using the full-waveform LiDAR by Equation (5) according to the parameters of the laser and the distances between the objects and the laser.
principle of the Figure 1, an overlapping area should exist between the two BFWSs detected by the APD 1 and the APD 2 to ensure that the proposed system works properly. Meanwhile, given that the complexity of the proposed system, inconsistent splitting ratio of the BS3, and the individual difference of the two APDs (APD1 and APD2), it is challenging to maintain the consistent characteristics of the two beams split by the BS 2. Therefore, the inconsistency of the two beams should be eliminated prior to utilization.

Differential Distance Selection
An overlapping area should exist between the two BFWSs detected by the two APDs. Therefore, the differential time (2L/c) should not exceed the minimum value of the received pulse width (τrmin). Thus, the differential distance should not exceed c/2 × τrmin, shown in Figure 7a. The differential BFWS of the two APDs is shown in Figure 7b when the differential distance is larger than c/2 × τrmin. The Figure 7b shows that the number of the NGZCPs is two, i.e., the number of the object is two. In fact, three objects are included in the simulated scene. Therefore, a differential distance that is larger than c/2 × τrmin may cause misjudgment of the number of the objects. Moreover, the differential BFWS is distorted, which increases the difficulty of wave decomposition and decreases the differential Gaussian fitting accuracy. In terms of the minimum value of the received pulse width τrmin, it can be evaluated before using the full-waveform LiDAR by Equation (5) according to the parameters of the laser and the distances between the objects and the laser.

Inconsistent Elimination of Two Beams
According to the principle of the differential optical path method in Figure 1, the BFWS is divided into two beams. The two beams are then detected by two different APDs. However, the amplitudes of the two BFWSs detected by the two APDs usually are different owing to the inconsistent splitting ratio of the BS 2 and individual difference of the two APDs. Under this situation, two variable gain amplifiers (VGAs) can be used to eliminate the inconsistency through gaining the BFWSs before entering the SC, as shown in Figure 8. The two amplitudes of the two BFWSs can be adjusted to the same desired value using the two VGAs, and the amplification coefficients are the ratios of the desired value and the two amplitudes of the two BFWSs. It is worth noting that the amplitudes are amplified when the differential BFWS is decomposed to obtain the valuable parameters. Therefore, the amplification coefficients should be divided prior to utilization so that the actual physical properties of the objects can be obtained.

Inconsistent Elimination of Two Beams
According to the principle of the differential optical path method in Figure 1, the BFWS is divided into two beams. The two beams are then detected by two different APDs. However, the amplitudes of the two BFWSs detected by the two APDs usually are different owing to the inconsistent splitting ratio of the BS 2 and individual difference of the two APDs. Under this situation, two variable gain amplifiers (VGAs) can be used to eliminate the inconsistency through gaining the BFWSs before entering the SC, as shown in Figure 8. The two amplitudes of the two BFWSs can be adjusted to the same desired value using the two VGAs, and the amplification coefficients are the ratios of the desired value and the two amplitudes of the two BFWSs. It is worth noting that the amplitudes are amplified when the differential BFWS is decomposed to obtain the valuable parameters. Therefore, the amplification coefficients should be divided prior to utilization so that the actual physical properties of the objects can be obtained.

Conclusions
A novel de-noising method for full-waveform LiDAR based on the differential optical path is proposed, and the mathematical models of the proposed method are developed and verified. On the basis of the proposed method, simulations are carried out, including the BGN, stop timing moment discrimination, the SNR, the waveform decomposition, and the differential Gaussian fitting. The following conclusions are obtained. (1) The proposed method can effectively suppress the BGN. (2) The BFWS is detected by two APDs placed before and after the focus of the focusing lens. Hence, the PP is transformed into NGZCP, which is more beneficial for the stop timing moment discrimination for the timer. (3) The SNR of the proposed method is improved but the improvement gets smaller by increasing distance. The relative increment percentage ΔSNR of the proposed method increases dramatically from 65 to 145 when the distance increases from 300 m to 1000 m, and it increases slightly from 145 to 161 when the distance increases from 1000 m to 3000. (4) The differential Gaussian fitting based on the LM algorithm can decompose the differential BFWS with high accuracy. The maximum absolute errors of the three objects are 3.56 × 10 −12 , 2.80 × 10 −11 , and 2.08 × 10 −12 W, respectively. The minimum absolute errors of the three objects are −6.23 × 10 −12 , −1.99 × 10 −11 , and −4.58 × 10 −12 W, respectively. (5) The differential distance should not be larger than c/2 × τrmin, and employing two VGAs can eliminate inconsistency of the two beams. The proposed method is applicable for the fullwaveform LiDAR application fields, such as surface topography and airborne vegetation mapping. In this study, our works mainly focus on the theoretical framework and simulation experiments for validating the proposed method. In our future works, we will carry out the experiments and verify the validation of the proposed method.

Conclusions
A novel de-noising method for full-waveform LiDAR based on the differential optical path is proposed, and the mathematical models of the proposed method are developed and verified. On the basis of the proposed method, simulations are carried out, including the BGN, stop timing moment discrimination, the SNR, the waveform decomposition, and the differential Gaussian fitting. The following conclusions are obtained. (1) The proposed method can effectively suppress the BGN.
(2) The BFWS is detected by two APDs placed before and after the focus of the focusing lens. Hence, the PP is transformed into NGZCP, which is more beneficial for the stop timing moment discrimination for the timer. (3) The SNR of the proposed method is improved but the improvement gets smaller by increasing distance. The relative increment percentage ∆SNR of the proposed method increases dramatically from 65 to 145 when the distance increases from 300 m to 1000 m, and it increases slightly from 145 to 161 when the distance increases from 1000 m to 3000. (4) The differential Gaussian fitting based on the LM algorithm can decompose the differential BFWS with high accuracy. The maximum absolute errors of the three objects are 3.56 × 10 −12 , 2.80 × 10 −11 , and 2.08 × 10 −12 W, respectively. The minimum absolute errors of the three objects are −6.23 × 10 −12 , −1.99 × 10 −11 , and −4.58 × 10 −12 W, respectively. (5) The differential distance should not be larger than c/2 × τ rmin , and employing two VGAs can eliminate inconsistency of the two beams. The proposed method is applicable for the full-waveform LiDAR application fields, such as surface topography and airborne vegetation mapping. In this study, our works mainly focus on the theoretical framework and simulation experiments for validating the proposed method. In our future works, we will carry out the experiments and verify the validation of the proposed method.