Point cloud noise reduction method and smart noise reduction system for photon counting lidar

Abstract

This invention provides a point cloud noise reduction method and a smart noise reduction system for photon counting lidar with adaptive optimization based on geographic semantic information. [Solution] The method synchronously acquires an optical image of the target region and the initial photon point cloud. After identifying the dominant feature type in the laser spot region by real-time semantic partitioning, the scene complexity is calculated based on the point cloud data and its grade is determined. Subsequently, based on the combination of feature type and scene complexity grade, a pre-set 2D parameter mapping table is queried to dynamically acquire the optimal denoising parameter set. Finally, adaptive denoising processing is driven using this parameter set to calculate the local density and directional consistency of points and construct a 2D discrimination space, thereby realizing the identification of signal points and noise points.

Description

[Technical Field]

[0001] The present invention belongs to the technical field of photon lidar data processing, and more specifically relates to a point cloud noise reduction method and smart noise reduction system for photon counting lidar with adaptive optimization based on feature semantic information. [Background technology]

[0002] Photon-counting LiDAR (Photon-Counting LiDAR) is an advanced active optical telemetry technology that detects single-photon echo signals, significantly reducing system power consumption and volume, enabling lightweight, compact, and highly sensitive detection. It is particularly well-suited for high-speed mobile surveying when mounted on unmanned aerial vehicles. However, its extremely high detection sensitivity has both advantages and disadvantages; while it receives weak laser echo signals, it inevitably introduces a large amount of background noise. This noise appears as discrete and randomly distributed noise points in the point cloud data, and its density and distribution are closely related to factors such as feature type, ambient lighting, and water turbidity, resulting in significant differences in the signal-to-noise ratio across different regions.

[0003] Chinese Patent Application Publication No. 117572381 discloses a single-photon lidar data processing method based on spatiotemporal correlation. While this invention achieves noise reduction and improved signal set and enhances the utilization of echo photons, this purely software method processes the point cloud only from a mathematical and statistical standpoint and therefore cannot sense or identify physical sources of noise (e.g., sudden changes in sunlight, platform vibrations, changes in water turbidity, etc.), and thus cannot employ targeted suppression strategies for specific types of noise.

[0004] Chinese Patent Application Publication No. 119335556 discloses a technology that enables real-time elevation imaging and high-precision imaging of underwater scenes using an underwater single-photon SPAD array imaging system. However, this system cannot perform real-time positioning, making it difficult to exert a significant impact in the field of topographic surveying.

[0005] The paper "Luan KF, Zhang KN, Qiu ZG, Wang J, Wang ZH, Xue Y, Zhu WD, Ling DD and Zhao X Y. 2023. Denoising method for light weight photon counting LiDAR based on an improved local sparse coefficient. National Remote Sensing Bulletin, 27(2):520-532" first calculates the effective signal interval of the initial photon observation data using a histogram-based statistical method, then applies coarse denoising to the data within the interval using a grid statistical method, and finally achieves fine denoising of the initial photon observation data using an improved local sparse coefficient method. Because this algorithm is complex in order to achieve high-precision denoising, the computational cost is enormous, making it difficult to implement real-time processing on embedded platforms with limited computational resources, such as unmanned aerial vehicles, thus limiting surveying efficiency.

[0006] As described above, existing technologies have the following main limitations: fixed parameters and lack of scene adaptability; a single discrimination dimension making it difficult to solve the problem of confusing signals and noise; a significant contradiction between algorithmic complexity and real-time capabilities; and poor physical sensing ability. Therefore, there is a strong need for a completely new point cloud denoising method for photon counting lidar devices that intelligently senses the scene, adaptively sets parameters, efficiently distinguishes signals and noise, and meets the requirements of real-time processing. [Prior art documents] [Patent Documents]

[0007] [Patent Document 1] Chinese Patent Application Publication No. 117572381 Specification [Patent Document 2] Chinese Patent Application Publication No. 119335556 Specification [Overview of the project] [Problems that the invention aims to solve]

[0008] The objective of the present invention is to solve the problem that conventional photon counting lidars process all scenes with the same set of fixed parameters, resulting in serious signal erroneous deletion in sparse scenes such as vegetation, and significant noise retention in dense scenes such as buildings, making it impossible to achieve balanced and optimized noise reduction performance in cross-feature scenes. [Means for solving the problem]

[0009] To achieve the above objective, the present invention provides a method for removing point cloud noise in a photon counting lidar, comprising the following steps. Synchronized acquisition step: This step includes synchronously acquiring an optical image of the target region and the original photon point cloud data collected by a photon counting lidar.

[0010] Feature type identification step: Perform real-time semantic partitioning on the optical image to identify the dominant feature type in the region covered by the laser spot in the original photon point cloud data.

[0011] Scene complexity determination step: Based on the original photon point cloud data, a statistic is calculated that includes at least one of the following indicators: point density variance, elevation relief, and noise ratio estimate. Based on this statistic, the scene complexity grade of the current region is determined.

[0012] Adaptive parameter acquisition step: Based on the dominant feature type and scene complexity rating, a two-dimensional parameter mapping table stored in a memory device is queried, and a denoising parameter set corresponding to the "feature type-scene complexity rating" combination is dynamically acquired. Here, the two-dimensional parameter mapping table stores optimized denoising parameters for each combination of different feature types and different scene complexity ratings.

[0013] Adaptive denoising execution step: Adaptive denoising is performed on the original photon point cloud data using the denoising parameter set. Here, the adaptive denoising process includes the following: For each point in the point cloud, calculate its local density ρ and the density gradient direction consistency C with neighboring points; construct a two-dimensional discriminant space based on the local density ρ and direction consistency C; and discriminate each point as a signal point or a noise point according to its position in the two-dimensional discriminant space.

[0014] In a preferred embodiment, when the adaptive noise reduction execution step determines whether a signal point or a noise point based on the position of each point in a two-dimensional discrimination space, the following determination rules are specifically included: The local density ρ of the point is equal to the first density threshold ρ min The above conditions are met, and the directional consistency C is equal to the consistency threshold C. min If the above conditions are met, the point is determined to be a core signal point and retained. The local density ρ of the point is equal to or greater than the second density threshold, and its directional consistency C is equal to the consistency threshold C. min The above conditions are met, and the gradient amplitude is the gradient threshold g max If it is greater than the first density threshold ρ, the point is determined to be an edge signal point and is retained. Here, the second density threshold is equal to the first density threshold ρ min It is assumed to be smaller. The local density ρ of a point is the first density threshold ρ min The above is true, but the directional consistency C is the consistency threshold C minIf the value is less than the second density threshold, the point is determined to be a high-density noise cluster and is removed. If the local density ρ of a point is less than the second density threshold, the point is determined to be low-density background noise and is removed. This embodiment solves the problem of difficulty in distinguishing between high-density noise clusters and true edge signals. Conventional density clustering methods rely on only a single density threshold and cannot effectively distinguish between high-density but directionally congested noise clusters formed by solar flares, specular reflections, etc., and true edge signals that are located at feature boundaries, are low-density, but have a consistent direction, leading to misidentification of the two.

[0015] In a preferred embodiment, the adaptive denoising execution step specifically includes the following substeps, which are performed sequentially: Coarse denoising substep: Preprocessing is performed on the original photon point cloud data using a fast filtering method based on statistical model fitting to remove discrete outlier noise points and obtain an initial filtered point cloud. Fine denoising substep: Adaptive denoising of a two-dimensional discriminant space based on local density ρ and directional consistency C is performed on the initial filtered point cloud. This embodiment significantly reduces the amount of input data for the fine denoising algorithm and dramatically improves the overall processing efficiency, effectively solving the technical problem that existing high-precision denoising algorithms have high computational costs and are difficult to meet the requirements for real-time processing on embedded platforms with limited computational resources such as unmanned aerial vehicles.

[0016] In a preferred embodiment, a fast filtering method based on statistical model fitting is specifically a statistical filtering method based on single Gaussian model fitting and a 3σ criterion, and includes the following steps: Divide the original photon point cloud data into multiple time fragments along the direction of the lidar's track; perform a single Gaussian model fitting on the set of elevation values ​​of the photon data within each time fragment, and calculate the mean μ and standard deviation σ of the photon elevation values ​​within that fragment, based on a 3σ criterion, retain only the photon data whose elevation values ​​fall within the interval [μ-2σ, μ+2σ], and constitute an initial filtered point cloud.

[0017] In a preferred embodiment, before or simultaneously with the synchronization acquisition step, a hardware pre-filtering step is further included. Based on the approximate distance information of the target area acquired by the auxiliary laser ranging module, the central processing module dynamically calculates a time gating window, and uses the time gating window to control the main detection unit of the photon counting lidar, so that photon data is collected only in the time period of the predicted target echo, thereby realizing pre-filtering at the hardware level for background noise. This step suppresses background noise at non-target distances from the physical aspect, and realizes targeted suppression of the physical source of noise. It complementarily solves the technical problem that pure software algorithms only start from the statistical characteristics of point cloud data, cannot perform source suppression based on the physical mechanism of noise generation (for example, background light noise at different distances), and thus there is a bottleneck in noise removal capability.

[0018] In a preferred embodiment, in the scene complexity determination step, the statistic is a scene complexity index, which is derived by fusing quantized values of three indicators: point cloud density variance, elevation relief and estimated noise ratio. Based on the value of the scene complexity index, the current area is determined as one of three predefined scene complexity grades: simple, medium or complex. The scene complexity index SCI is calculated by the formula: SCI=w1·normalize(σ density )+w2·normalize(H variance )+w3·R noise , where σ density is the variance value of point cloud density, H variance is the elevation relief, and R noiseis the noise ratio estimate, w1, w2, and w3 are weighting coefficients, and normalize(·) represents the operation of normalizing to the [0,1] interval. Simple scenes, moderate scenes, and complex scenes correspond to cases where the Scene Complexity Index SCI is less than the first threshold, between the first and second thresholds, and greater than or equal to the second threshold, respectively. Elevation relief H variance Formula: H variance =α·normalize(H range )+β·normalize(σ height This is calculated by ), where H range σ is the difference in elevation of the target region. height is the standard deviation of elevation, α and β are weighting coefficients, and normalize(·) represents the normalization operation.

[0019] In a preferred embodiment, the two-dimensional parameter mapping table queried in the adaptive parameter acquisition step includes, for each combination of feature type and scene complexity grade, the following parameters: number of neighbors k, horizontal search scale a and vertical search scale b of the anisotropic Gaussian kernel, Gaussian kernel bandwidth parameter h, and minimum density threshold ρ. min , minimum direction consistency threshold C min , maximum gradient amplitude threshold g max The key feature is that at least one of the following differential values ​​is set: the threshold for the number of points near the edge point core, and the tolerance for the pulse elevation error ΔH.

[0020] The present invention provides a smart noise reduction system for realizing the point cloud noise reduction method, which includes: an optical imaging module for performing optical image acquisition in a synchronous acquisition step; a photon counting lidar main detection unit for performing initial photon point cloud data acquisition in a synchronous acquisition step; and a central processing module connected to the optical imaging module and the main detection unit, configured to perform a feature type identification step, a scene complexity determination step, an adaptive parameter acquisition step, and an adaptive noise reduction execution step. [Effects of the Invention]

[0021] 1. This invention fundamentally solves the performance imbalance problem of fixed-parameter methods through a two-dimensional adaptive mechanism of "feature type-scene complexity" that achieves cross-scene balanced high-precision noise reduction, maintaining excellent noise reduction performance in both vegetation and building areas, and significantly improving the fidelity of topographic features.

[0022] 2. This invention introduces density gradient directional consistency (C) to overcome the challenge of distinguishing between high-density noise and edge signals, and constructs a two-dimensional discrimination space together with density (ρ). This fundamentally distinguishes between high-density noise clusters and sparse edge signals, significantly reducing the misclassification rate.

[0023] 3. The present invention employs a collaborative architecture of "hardware pre-filtering + software hierarchical processing (coarse noise reduction + fine noise reduction)" to meet the real-time processing requirements of embedded platforms, significantly reducing the computational load and enabling the high-precision noise reduction algorithm to operate in real time on resource-constrained platforms such as unmanned aerial vehicles.

[0024] 4. The present invention constructs an intelligent and fully automated processing flow, eliminating the need for human intervention throughout the entire process, from semantic sensing of optical images to point cloud complexity analysis, dynamic parameter setting, and algorithm execution, allowing the system to autonomously adapt to complex and variable real-world surveying environments. Other advantages, purposes, and features of the present invention will be partially revealed by the following description and partially understood by those skilled in the art through the study and practice of the present invention. [Brief explanation of the drawing]

[0025] [Figure 1] This is an overall flowchart of the point cloud noise reduction method in a photon counting lidar according to the present invention. [Figure 2] This is a schematic diagram illustrating the principle of semantic identification of features and adaptive parameter setting according to the present invention. [Figure 3] This is a logic module configuration diagram of the smart noise reduction system according to the present invention. [Figure 4] This figure shows the results of a comparative experiment between the noise reduction method of the present invention and the DBSCAN method using fixed parameters. [Modes for carrying out the invention]

[0026] The present invention will be described in more detail below with reference to examples, so that those skilled in the art can implement it based on the description in the specification.

[0027] The photon counting lidar point cloud noise reduction method according to the present invention is an adaptively optimized photon counting lidar point cloud noise reduction method based on feature semantic information. First, refer to the overall flowchart of the present invention shown in Figure 1 and the adaptive closed-loop principle shown in Figure 2. By sequentially executing the following steps, a complete sense-decision-execute intelligent closed-loop is formed, which specifically includes the following steps.

[0028] Step S1: Synchronization acquisition step The system controls a synchronized optical camera fixed to the optical axis of the photon-counting lidar. At the same precise time that the main detection unit of the photon-counting lidar emits a single laser detection pulse, the optical camera is triggered to expose the same target area and collect a single frame of optical image. Simultaneously, the main detection unit receives the laser echo and generates primal photon point cloud data that is spatially precisely corresponding to the optical image and temporally synchronized, via Time-Correlated Single Photon Counting (TCSPC) technology.

[0029] Step S2: Feature type identification step The optical image collected in step S1 is input into a pre-trained semantic segmentation model for real-time inference. This model identifies the feature category to which each pixel in the image belongs. In this invention, we focus on the dominant land cover type of the area covered by the current laser spot, such as open bodies of water, dense vegetation (e.g., forests, shrubs), artificial structures, bare ground, etc.

[0030] In the specific implementation of the present invention, the semantic partitioning model preferentially employs a U-Net convolutional neural network based on an encoder-decoder structure. Of these, the encoder (i.e., the downsampling path) typically includes 4 to 5 stages, each stage consisting of a convolutional layer, a ReLU activation function, and a max pooling layer, used to extract and compress features stepwise, and the number of output channels increases stepwise with the progression of the stages. The decoder (i.e., the upsampling path) has a structure symmetric to the encoder, achieving upsampling through transposed convolution or interpolation operations, fusing the feature maps of the corresponding stages of the encoder using "skip connections," and further combining deep semantic information with shallow detailed information. The end of the network is a 1x1 convolutional layer, the number of output channels of which is equal to the total number of feature categories C, and a Softmax function is connected after this convolutional layer to finally obtain the category probability distribution of each pixel.

[0031] An annotated aerial imagery dataset containing diverse typical survey scenes is used to train the semantic segmentation model. To resolve the category pixel imbalance problem, a weighted cross-entropy loss function is selected during the training process. The average loss value L of the batch data is calculated according to JPEG0007914992000002.jpg689. The definitions of each parameter in the formula are as follows: L represents the loss value, C represents the total number of feature categories, and y c This is an indicator variable (value is 0 or 1) used to indicate whether the true category of the current pixel is c. JPEG0007914992000003.jpg6150ω c This is a weight coefficient specifically set for category c. Weight ω c The learning rate is calculated inversely based on the frequency of occurrence of each category in the training dataset, thereby enhancing the model's learning effect on rare categories. The training process unfolds via a backpropagation algorithm, employing the Adam optimizer and iteratively optimizing the model's weights with typical parameters such as an initial learning rate of 1e-4.

[0032] In this invention, the input to the semantic segmentation model is an optical image collected in strict synchronization with the photon point cloud, thereby ensuring high-precision alignment of the two in the spatiotemporal dimensions. After the model outputs a pixel-level category map, the dominant feature type (e.g., dense vegetation) of the region corresponding to the laser spot is statistically determined, and this dominant feature type is used as one of the core inputs to query a "feature type-scene complexity" 2D parameter mapping table, which then drives a subsequent adaptive denoising flow, achieving a tight coupling between semantic sensing and point cloud processing.

[0033] It should be noted that the original features of this invention do not lie in the semantic segmentation model or its learning method itself. These are known technical means for realizing the function of "feature type identification." The core inventive aspect of this invention lies in the fact that it utilizes the results of the above feature type identification to refer to a "feature type-noise reduction parameter mapping table" and realizes adaptive parameter setting, in which the completely novel method and control logic of this invention reside.

[0034] Step S3: Scene complexity determination step and adaptive parameter acquisition step This step is the core of achieving intelligent adaptability in this method, and its innovation lies in constructing a dual-axis integrated decision mechanism based on "feature type" and "scene complexity," replacing conventional single-feature mapping or fixed-parameter models. Specifically, it includes the following:

[0035] S3.1 Quantitative Evaluation of Scene Complexity After identifying the feature types, it is necessary to further quantify the complexity of the current point cloud scene. The following quantitative analysis process is performed on the raw photon point cloud data acquired in step S1. 1. Divide the current area covered by the laser spot into multiple regular spatial grids (for example, a 10m x 10m grid). 2. Perform feature extraction, specifically including the following: The variance of the point cloud density (σ density First, the average point density within each grid is calculated, and then the variance of the point densities across all grids is statistically analyzed. This value reflects the heterogeneity in the spatial distribution of the point cloud. Elevation relief (H variance ): In order to comprehensively represent the topographic relief, formula H variance =α·normalize(H range )+β·normalize(σ height ) is used to calculate it. Here, H range σ represents the difference in elevation within the current region (i.e., the difference between the maximum and minimum elevations), height The value represents the standard deviation of the elevation values ​​for all points. α and β are weighting coefficients (for example, α=0.6 and β=0.4; however, the present invention is not limited to these values), and normalize(·) represents the operation of normalizing the values ​​to the interval [0,1]. Noise ratio estimate (R noise ): To quickly evaluate the noise level, a certain number of sample points (e.g., 1000) are randomly selected from the current point cloud. For each sample point, the average distance to its k nearest neighbors is calculated, and then the proportion of all sample points whose average distance exceeds a pre-set empirical threshold is statistically calculated and this is expressed in R noise This value reflects the proportion of discrete and anomalous noise points. Calculation and Classification of Scene Complexity Index (SCI): The three feature indices described above are weighted and merged to obtain the Scene Complexity Index (SCI). The Scene Complexity Index SCI is given by the formula: SCI = w1·normalize(σ density )+w2·normalize(H variance )+w3·R noise This is expressed as follows: Here, w1, w2, and w3 are pre-set weighting coefficients (e.g., w1=0.35, w2=0.30, w3=0.35). Based on the SCI value, the scene complexity is classified into three predefined grades: simple scenes (SCI<0.35, e.g., open, flat bodies of water, uniform bare ground), moderate scenes (0.35≦SCI<0.70, e.g., suburban vegetation, low-density building areas), and complex scenes (SCI≧0.70, e.g., densely populated urban areas, mountainous vegetation, land-water boundary zones).

[0036] S3.2 Dynamic parameter acquisition based on a 2D mapping table The system of the present invention has a pre-configured "feature type-scene complexity" 2D parameter mapping table (see example in Table 1 below). This mapping table is built through a large number of preliminary experiments and stores optimized denoising parameters for different "feature type-scene complexity grade" combinations. Based on the feature type identified in step S2 and the scene complexity grade determined in step S3.1, the system automatically queries the 2D parameter mapping table and dynamically obtains the corresponding optimal denoising parameter set.

[0037] [Table 1]

[0038] Explanation of the parameters in Table 1 Number of neighbors in k:k nearest neighbor search a,b: Horizontal and vertical search scales (meters) of anisotropic Gaussian kernels. h: Gaussian kernel bandwidth parameter ρmin : Minimum density threshold C min : Minimum Directional Consistency Threshold g max : Maximum gradient amplitude threshold Edge Neighborhood Count: The minimum number of core neighbors required to hold an edge point. ΔH: Permissible error in the same pulse elevation (meters)

[0039] Step S4: Adaptive noise reduction execution step In step S3.2, the dynamically acquired parameter set is injected into the denoising algorithm, denoising is performed on the original photon point cloud, and finally, a high-precision point cloud after denoising is output. In this invention, in order to achieve both processing efficiency and denoising accuracy, a hierarchical denoising strategy from coarse to fine is given priority, and specifically includes the following.

[0040] S4.1 Coarse Noise Reduction – Statistical-Based High-Speed ​​Filtering To improve overall processing efficiency, we first perform high-speed coarse filtering on the initial photon point cloud data, primarily to remove obvious outlier noise points. The specific procedure is as follows: Time interval division: This involves dividing the point cloud data into continuous time intervals along the direction of the LiDAR track (for example, each time interval corresponds to 0.1 seconds of data). Fitting a single Gaussian model: Set of elevation values ​​for all photons within each time interval H={h1,h2,...,h n Statistical analysis is performed on the data to calculate the mean μ and standard deviation σ of the photon elevation values ​​within the given time interval. The calculation formulas are as follows: μ = (1 / n) × Σh i σ = √[(1 / n) × Σ(h i -μ) 2 ]. 3σ-based filtering: Based on the assumption of a normal distribution, a 3σ criterion is applied, retaining only photon data where the elevation value falls within the interval [μ-2σ, μ+2σ] (theoretically, this interval contains approximately 95.4% of normally distributed data). This step allows for the rapid removal of a large number of discrete noise points that clearly deviate from the elevation distribution of the main subject, thereby obtaining an initial filtered point cloud and significantly reducing the computational load of subsequent refinement.

[0041] S4.2 High-precision noise reduction—Adaptive density-gradient fusion clustering (ADGFC) The Adaptive Density-Gradient Fusion Clustering Algorithm (ADGFC), which is the core of the present invention, is executed on the initial filtered point cloud described above. The core innovation point of this algorithm is the construction of a two-dimensional discrimination space of "density (ρ) - directional consistency (C)," which enables high-precision differentiation of signals and noise that were difficult to distinguish using conventional methods.

[0042] Any point P in the point cloud i (x i ,y i ,z i For ), perform the following steps: Step 1: Adaptive neighborhood density calculation First, based on the k value set for the scene in the parameter set, point P is created using a spatial index structure (e.g., KD-Tree). i We quickly search for the k nearest neighbors and construct the neighborhood set N_k(i). Next, P i and each neighboring point P j To quantify the spatial proximity to [the object], an anisotropic Gaussian kernel function is used to calculate the weights, and the formula is as follows: W(P i ,P j )=exp(-D_aniso 2 (P i ,P j ) / (2h 2 )) Here, the core part is the anisotropy distance squared D_aniso 2 The formula is: D_aniso 2 =(Δx 2 +Δy 2 ) / a 2 +Δz 2 / b 2 It is defined as follows. The coordinate difference between two points is Δx = x i -x j Δy=yi -y j , Δz=z i -z j , wherein a and b are the horizontal search scale and vertical search scale in the adaptive parameter set respectively, and h is a bandwidth parameter obtained by referring to the feature type. By introducing these two independent scale parameters, the distance measurement can be flexibly adapted to the spatial structure characteristics of different features. That is, for a building surface extending in the horizontal direction, a small a value converges the search range in the horizontal direction; for vegetation growing in the vertical direction, a small b value converges the search range in the vertical direction, so that the definition of neighborhood is more consistent with the true geometric morphology of the feature. h is the bandwidth parameter of the Gaussian kernel, which is also specified by the parameter set. Finally, point P i has local density ρ(P i ) which is calculated by summing the Gaussian weights of all its neighboring points, and the formula is as follows: ρ(P i )=Σ_{j∈N k (i)}w(P i ,P j ). This numerical value reflects the compactness with which point P i is surrounded by other points in the local space defined by anisotropy (where N k (i) represents the k-nearest neighbor set of P i ).

[0043] Step 2: Construction of density gradient field After obtaining the local density, the algorithm further constructs a density gradient field to analyze the spatial change trend of density. First, calculate the density gradient vector of point P i . The formula is: ∇ρ(P i )=Σ_{j∈N k}(i)}[(P j -P i )·(ρ(P j )-ρ(P i ))·w(P i ,P j )] / ||P j -Pi|| 2 . This vector is a vector pointing to each neighboring point (P j -P i) is obtained by weighted composition, and the weight is the density difference ρ(P j )-ρ(P i )) and Gaussian weight w(P i ,P j ) is directly proportional to. ∇ρ(P i ) indicates the spatial direction in which the local density at point P i increases most rapidly. For points on the surface of a real feature, the gradient direction usually points to the density center inside the feature body, while for random noise points, the gradient direction is generally disordered. Next, calculate the amplitude (i.e., norm) of the gradient vector. The formula is:g(P i )=||∇ρ(P i )||. This is an index representing the intensity of density change. Finally, when g(P i )>0, calculate the unit vector in the normalized gradient direction. The formula is:n(P i )=∇ρ(P i ) / g(P i ). This is used for the subsequent direction consistency analysis.

[0044] Step 3: Measurement of Direction Consistency To quantify the coordination between a certain point and its surrounding environment, the algorithm performs the measurement of direction consistency. First, for point P i and any of its neighboring points P j , calculate the inner product of their normalized gradient direction vectors, and take this as the direction consistency c(P i ,P j ) between the two points. The formula for c(P i ,P j ) is: c(P i ,P j )=n(P i )·n(P j ). The result of the inner product is in the range of [-1,1]. The closer the result is to 1, the more consistent the density increase directions at the two points are, and the higher the probability that they belong to the same continuous feature. The closer the result is to 0, the more approximately orthogonal the directions are, and the point is more likely to be located at the feature boundary of the feature. The closer the result is to -1, the opposite the directions are, and the higher the probability that one of the points is in a noise region. Since gradients in opposite directions are often caused by noise interference, the average direction consistency C(P of point P i ​i When calculating C(P), we ignore the negative contribution, and the formula is as follows: i )=(1 / k)·Σ_{j∈N_k(i)}max(0,c(P i ,P j )) C(P i ) is a scalar in the range of 0 to 1, and the higher its value, the higher the point P i The direction of density changes in the neighboring region is more unified and coordinated, meaning that the probability of that point belonging to a truly ordered structure is high.

[0045] Step 4: Classification of point types based on a 2D discriminant space The algorithm utilizes three key scalars (local density ρ, mean directional consistency C, and gradient amplitude g) calculated from a parameter set and dynamically obtains a threshold (minimum density threshold ρ) from the parameter set. min , minimum direction consistency threshold C min , maximum gradient amplitude threshold g max ) combined with the point P in the (ρ,C) 2D discriminant space i It classifies them with high accuracy. The rules are as follows: 1. Core signal point: ρ(P i )≧ρ min and C(P i )≧C min If P i These are identified as high-density, high-matching core signal points and directly retained. 2. Edge signal point: ρ(P i )≧0.6·ρ min , C(P i )≧C min and g(P i )>g max If P i These points are identified as potential edge signal points. Although the density of these points does not meet the criteria for core signal points, their direction is coordinated with the neighbors, and they are located in areas with rapid density changes (high gradients), which matches the characteristics of feature edges, so we proceed to the next step, secondary verification. 3. High-density noise population: ρ(P i )≧ρ min But C(P i ) <C minIf P i These are identified as high-density noise points. Although these points have a high local density, their directional consistency is extremely low and their direction is disordered. Therefore, they match the characteristics of cohesive noise such as bright spots in sunlight and false images of specular reflections, and are removed. 4. Low density background noise: ρ(P i )<0.6·ρ min If P i These are identified as sparse background noise points and deleted directly.

[0046] Step 5: Secondary verification of edge points To minimize the erroneous removal of weak signals from true feature boundaries for points classified as "edge signal points," the algorithm performs a secondary verification. The specific operation is as follows: The number of points already identified as "core signal points" within the k-neighborhood set of the edge point is statistically calculated, and this is then calculated as N core This value is then compared to the "number of edge neighbors" threshold in the parameter set. core If the threshold is exceeded, it means that there are sufficient core signal points acting as "support" at that edge point, and it is determined to be part of a true feature and ultimately retained. Otherwise, it is determined that there is insufficient support and that it is likely isolated noise, so it is removed. This step significantly improves signal retention and topographic integrity in complex edge regions.

[0047] Step 6: Time Dimension Constraints To further leverage the physical properties of the photon-counting lidar and improve the reliability of the discrimination, the algorithm optionally adds a time-dimension constraint. First, the entire point cloud is grouped by emission pulse based on the laser pulse ID recorded in each photon data. For all echo points belonging to the same laser pulse (i.e., one pulse group), the following process is performed: If there is at least one point determined to be a "core signal point" in the above step within that pulse group, it is determined that the laser pulse successfully detected a true target. In this case, all points within that pulse group whose elevation values ​​are within a reasonable range are retained. The "reasonable range" is defined by the pulse elevation tolerance ΔH in the parameter set. For example, a range of ΔH above and below the elevation value of the core signal point within the pulse group, or the highest and lowest elevation values ​​within the group, is allowed. ΔH is preset based on the physical properties of the features. For example, ΔH is set large for bodies of water (e.g., 15m) because the elevation of echoes is easily dispersed due to subsurface scattering, while ΔH is set small for building surfaces (e.g., 5m) because echoes are concentrated. This constraint takes advantage of the characteristic that echoes of true targets are temporally correlated, and effectively recovers weak echo points belonging to the same true target that might be mistakenly deleted by a single criterion due to their sparse spatial distribution, thereby enhancing the robustness of the algorithm.

[0048] To further improve system performance, a hardware-level pre-filtering step can be added before step S1, and the specific operation is as follows: Before emitting the main laser pulse, the auxiliary laser ranging module first emits a low-power ranging pulse to the same target area to quickly obtain the approximate distance value Rough of the target. The central processing module calculates the predicted echo time window based on Rough, and the predicted echo center time T rough = 2 × Rough / c (where c is the speed of light). Next, set the distance tolerance ΔR (e.g., ±5m) and calculate the corresponding time tolerance ΔT = 2 × ΔR / c. Finally, control the time-correlated single-photon counting (TCSPC) module and set the time window [Trough -ΔT,T rough Photon counts are collected only within the +ΔT range. This step suppresses a large amount of background noise (atmospheric scattering, reflections from non-target areas, etc.) that is inconsistent with the target distance from the data source. According to measured data, this step can reduce the initial data volume by more than 60%, thereby significantly reducing the load on subsequent data processing.

[0049] The present invention provides a smart noise reduction system for realizing the above method. Please refer to the system logic module configuration diagram shown in Figure 3. The system includes the following core modules. Optical imaging module: Used to perform step S1 and synchronously acquire optical images of the target region. Photon counting lidar main detection unit: Performs step S1, emits a laser and receives echoes, and is used to generate the initial photon point cloud data. Central Processing Module: Functions as the core of the system and incorporates the following functional units: Semantic Identification Unit: Loads a pre-trained semantic segmentation model and performs step S2 to identify the dominant feature type in the area covered by the laser spot. Adaptive mapping library: Stores a "feature type-scene complexity" 2D parameter mapping table, executes step S3, runs a query based on the input feature type and scene complexity grade, and outputs the corresponding optimal denoising parameter set. Noise reduction control unit: Receives the parameter set output from the adaptive mapping library, drives and executes the hierarchical noise reduction algorithm (coarse noise reduction and ADGFC fine noise reduction) in step S4, and completes the entire flow processing. In addition, the system can optionally be equipped with an auxiliary laser ranging module, which is used to implement the hardware pre-filtering function described in the preceding paragraph.

[0050] Comparative experiment To objectively and quantitatively verify the technical effects of the present invention, rigorous comparative experiments were conducted using an actual drone-mounted surveying platform. The experimental site was selected as a typical complex scene in a certain park. This area simultaneously includes a densely vegetated area with a complex structure and a land-based building area with distinct features, allowing for thorough verification of the adaptability and performance balance of the noise reduction algorithm to different types of features.

[0051] In the experiment, data was collected using the same drone-mounted photon counting lidar system, and the system parameters were set as follows: laser wavelength 532nm, pulse frequency 2kHz, flight altitude 10m, and flight speed 5m / s. For the control group, a DBSCAN density clustering noise reduction method with fixed parameters, widely used in the industry, was employed. These parameters (density threshold = 5, search radius = 1.0m) were "optimal fixed values" optimized through extensive prior experiments for land-based building scenes. In the experimental group, the adaptive smart noise reduction method based on feature semantics proposed by this invention was fully adopted.

[0052] The experimental results confirmed that in densely vegetated areas, the signals and noise from the original photon point cloud were heavily mixed (see Figure 4(a)). The processing results of the DBSCAN method with fixed parameters (see Figure 4(b)) revealed its inherent flaws. Specifically, the density threshold (=5) optimized for the building area was excessively strict for the naturally vegetation point cloud, which has a relatively sparse distribution. As a result, there was a shortage of points in the local vicinity, causing a large number of true vegetation echo signals to be misidentified as noise and removed. Ultimately, this resulted in clear "voids" and hydrophobicity in the point cloud, causing serious distortion of the structural characteristics of the vegetation.

[0053] Furthermore, because this method uses only a single density threshold as its discrimination criterion, it cannot identify noise clusters (e.g., cluster-like noise) where the local density has just reached the threshold but is essentially directionally disordered, resulting in the incorrect retention of noise in these areas.

[0054] Quantitative statistical data showed that the method could only remove about 60% of the noise in the area, leaving 40% of noise points remaining. The signal retention rate (SRR) was at least 58.6%, and the topographic accuracy (TF) was only 0.73.

[0055] In contrast, the method of the present invention demonstrated superior intelligent adaptive capabilities in the area. The system performed real-time semantic segmentation on synchronously acquired optical images, automatically identifying the area as "dense vegetation," and then, after real-time analysis of point cloud data, calculated a scene complexity index SCI=0.72, determining it to be a "complex scene." Based on the combination condition of "dense vegetation + complex scene," the system dynamically queried a pre-configured 2D parameter mapping table and loaded the corresponding optimal parameter set (e.g., {k=12, a=0.8m, b=0.5m, ρ min = 5.5, C min =0.70, g max =0.30, ...}). Subsequently, the Adaptive Density-Gradient Fusion Clustering (ADGFC) algorithm is executed based on these parameters, and its core innovation point lies in the (ρ,C)2D discriminant mechanism. For example, local density ρ = 4.2(ρ min (=5.5 or less), but the directional consistency C=0.73(C min For points where the local density ρ = 6.8(ρ) is greater than 0.70 and located in a high-gradient region, the algorithm determines that the point is cooperative with the surrounding environment and classifies and retains it as an "edge signal point". Conversely, for points where the local density ρ = 6.8(ρ) min (=5.5 or higher) but directional consistency C=0.38(C minFor points where the value is less than 0.70, the algorithm determines that the direction of that point is disordered, classifies it as a "high-density noise cluster," and removes it. This discrimination mechanism fundamentally achieves a highly accurate distinction between weak signals and cohesive noise. As shown in Figure 4(c), the method of the present invention effectively suppresses noise while preserving the complete structure and detailed features of the vegetation to the greatest extent possible. According to quantitative statistical data, the noise retention rate in the area was significantly reduced to 15%, the signal retention rate improved remarkably to 91.2%, and the topographic integrity reached 0.92.

[0056] In areas with onshore structures, the echo signals from the building surfaces are strong and relatively densely distributed, so the DBSCAN method using fixed parameters (these parameters are optimized specifically for this type of scene) achieved good basic noise reduction (see Figure 4(e)). The specific performance indicators are as follows: signal retention rate 89.5%, noise suppression rate 92.1%, and terrain integrity 0.93. The method of the present invention automatically sets a parameter set corresponding to "buildings + medium-difficulty scenes" by semantic identification (e.g., {k=8, a=0.7m, b=0.3m, ρ min =5.0, C min =0.80, ...}), similarly excellent noise reduction results were achieved (see Figure 4(f)). Furthermore, overall and slight superiority was achieved in each performance index: signal retention rate 90.8%, noise suppression rate 93.6%, and terrain integrity 0.95. This superiority mainly stems from the direction consistency (C) discrimination criterion in the ADGFC algorithm. This criterion makes it possible to identify and remove anisotropic edge noise points, which were difficult to process with fixed-parameter methods, with high accuracy.

[0057] The conclusions of the integrated comparison are shown in the table below.

[0058] [Table 2]

[0059] Experimental conclusion The DBSCAN method using fixed parameters exhibited significant performance imbalances between different feature types, with a difference in signal retention between vegetation areas and land-based building areas reaching approximately 31 points. In contrast, the present invention employs a two-dimensional adaptive mechanism of "feature type-scene complexity" to achieve balanced optimization of performance across multiple scenes, reducing the difference in signal retention between the two feature areas to less than 1 point. Overall, the present invention improved the average signal retention by 17 points, the average noise suppression rate by 15.8 points, and the average topography integrity by 14.6 points. This experiment demonstrated the remarkable effectiveness of the present invention in solving the major technical challenges of "performance imbalance in fixed-parameter methods" and "confusion of high-density noise and edge signals." Furthermore, the hierarchical noise reduction strategy from coarse to fine employed by the present invention, when combined with an optional hardware pre-filtering function, effectively ensures the real-time processing capability of the algorithm in embedded platforms such as drones, realizing intelligence throughout the entire flow from data sensing and intelligent decision-making to high-precision execution.

[0060] As described above, preferred embodiments of the present invention have been explained, but no embodiment shall be limited to the embodiments of the present invention in any way. Any simple changes, equivalent variations, and modifications to the above embodiments based on the technical substance of the embodiments of the present invention shall still be included within the scope of the technical scheme of the embodiments of the present invention.

Claims

1. A synchronized acquisition step involves synchronously acquiring an optical image of the target region and the initial photon point cloud data collected by a photon counting lidar, A feature type identification step involves performing real-time semantic segmentation on an optical image to identify the dominant feature type in the region covered by the laser spot in the original photon point cloud data, and A scene complexity determination step involves calculating a statistical quantity that includes at least one of the following indicators based on the initial photon point cloud data: point density dispersion, elevation relief, and noise ratio estimate, and determining the scene complexity grade of the current region based on the said statistical quantity. Adaptive parameter acquisition step: Based on the dominant feature type and scene complexity grade, query a two-dimensional parameter mapping table stored in a storage means in advance, which stores optimized denoising parameters for each combination of different feature types and different scene complexity grades, and dynamically acquire a set of denoising parameters corresponding to the combination of feature type and scene complexity grade. A method for denoising a photon counting lidar point cloud data, comprising: an adaptive denoising execution step, which performs an adaptive denoising process on the original photon point cloud data using the aforementioned denoising parameter set, the process including: calculating the local density ρ and the density gradient direction consistency C with neighboring points for each point in the point cloud; constructing a two-dimensional discriminant space based on the local density ρ and the direction consistency C; and discriminating each point as a signal point or a noise point based on its position in the two-dimensional discriminant space.

2. In the adaptive noise reduction execution step, when determining whether a signal point or a noise point is determined based on the position of each point in the two-dimensional discrimination space, The local density ρ of a point is the first density threshold ρ ｍｉｎ The above conditions are met, and the directional consistency C is equal to the consistency threshold C. ｍｉｎ If the above conditions are met, the determination rule is to determine that the point is a core signal point and retain it, The local density ρ of a point is equal to the first density threshold ρ ｍｉｎ The second density threshold is smaller than or equal to the second density threshold, and its directional consistency C is equal to the consistency threshold C. ｍｉｎ The above conditions are met, and the gradient amplitude is the gradient threshold g. ｍａｘ If it is greater, the determination rule is to determine that the point is an edge signal point and retain it, The local density ρ of a point is the first density threshold ρ ｍｉｎ The above is true, but the directional consistency C is the consistency threshold C ｍｉｎ If the value is less than the specified value, the rule is to determine that the point is a high-density noise cluster and delete it. A determination rule is in which, if the local density ρ of a point is less than the second density threshold, the point is determined to be low-density background noise and deleted. The method for removing point cloud noise from a photon counting lidar according to claim 1, characterized by including the following:

3. The adaptive noise reduction execution steps are performed sequentially. The initial photon point cloud data is preprocessed using a fast filtering method based on statistical model fitting to remove discrete outlier noise points, and a coarse denoising substep is performed to obtain an initial filtered point cloud. A refined denoising substep is performed on the initial filtered point cloud, which involves performing adaptive denoising of a two-dimensional discriminant space based on local density ρ and directional consistency C. A method for removing point cloud noise from a photon counting lidar according to claim 1 or 2, characterized by including the following:

4. The fast filtering method based on statistical model fitting in the coarse noise reduction substep is a statistical filtering method based on single Gaussian model fitting and a 3σ criterion. The process involves dividing the initial photon point cloud data into multiple time fragments along the direction of the LiDAR's flight path, For each time fragment, a single Gaussian model fitting is performed on the set of elevation values ​​of the photon data, and the mean μ and standard deviation σ of the photon elevation values ​​within that fragment are calculated. Based on the 3σ criterion, the process involves retaining only photon data whose elevation values ​​fall within the interval [μ-2σ, μ+2σ], and constructing an initial filtered point cloud. The method for removing point cloud noise from a photon counting lidar according to claim 3, characterized by including the following:

5. The method for denoising a point cloud of a photon counting lidar according to claim 4, further comprising a hardware-level pre-filtering step before or simultaneously with the synchronous acquisition step, wherein the central processing module dynamically calculates a time gating window based on approximate distance information of the target area acquired by the auxiliary laser ranging module, uses the time gating window to control the main detection unit of the photon counting lidar, and causes the main detection unit to collect photon data only during the predicted time period of the target echo, thereby performing hardware-level pre-filtering for background noise.

6. The method for removing point cloud noise from a photon counting lidar according to claim 1, characterized in that, in the scene complexity determination step, the statistical quantity is a scene complexity index, which is derived by fusing the quantitative values ​​of three indicators: point cloud density variance, elevation relief, and noise ratio estimate, and based on the numerical value of the scene complexity index, it is determined whether the current area is one of three predefined scene complexity grades, namely simple scene, moderate scene, or complex scene.

7. The scene complexity index SCI is calculated by the formula SCI=w1·normalize(σ ｄｅｎｓｉｔｙ )+w2·normalize(H ｖａｒｉａｎｃｅ )+w3·R ｎｏｉｓｅ , which σ ｄｅｎｓｉｔｙ is the variance value of the point cloud density, H ｖａｒｉａｎｃｅ is the elevation relief, R ｎｏｉｓｅ w1, w2, and w3 are the noise ratio estimates, and normalize(•) represents the operation of normalizing to the [0,1] interval. Simple scenes, moderate scenes, and complex scenes correspond to cases where the Scene Complexity Index (SCI) is below the first threshold, between the first and second thresholds, and above the second threshold, respectively. Elevation relief degree H ｖａｒｉａｎｃｅ This is the formula H ｖａｒｉａｎｃｅ =α・normalize(H ｒａｎｇｅ )+β・normalize(σ ｈｅｉｇｈｔ It is calculated by the following method: H ｒａｎｇｅ σ is the difference in elevation of the target region. ｈｅｉｇｈｔ The method for removing point cloud noise from a photon counting lidar according to claim 6, characterized in that is the standard deviation of elevation, and α and β are weighting coefficients.

8. In the adaptive parameter acquisition step, the 2D parameter mapping table queried provides, for each combination of feature type and scene complexity grade, the number of neighbors k, the horizontal search scale a and vertical search scale b of the anisotropic Gaussian kernel, the Gaussian kernel bandwidth parameter h, and the minimum density threshold ρ. ｍｉｎ , minimum direction consistency threshold C ｍｉｎ , maximum gradient amplitude threshold g ｍａｘ The method for removing point cloud noise from a photon counting lidar according to claim 1, characterized in that at least one differential value is set from among the parameters of the threshold value of the number of points near the edge point core and the pulse elevation tolerance ΔH.

9. The method for denoising a point cloud of a photon counting lidar according to claim 8, characterized in that when calculating the local density ρ in the adaptive denoising execution step, a weighted sum method based on an anisotropic Gaussian kernel function is adopted, and the anisotropic distance is determined based on the horizontal search scale a and the vertical search scale b in the denoising parameter set.

10. An optical imaging module for performing optical image acquisition in a synchronous acquisition step. A photon counting lidar main detection unit, which is for performing the collection of initial photon point cloud data in the synchronous acquisition step. An intelligent noise reduction system for realizing the point cloud noise reduction method according to claim 1, characterized in that a central processing module is connected to the optical imaging module and the main detection unit and is configured to perform a feature type identification step, a scene complexity determination step, an adaptive parameter acquisition step, and an adaptive noise reduction execution step.

Images