A laser radar ground point extraction system based on visual cognitive simulation filtering

Abstract

This invention discloses a LiDAR ground point extraction system based on visual cognitive simulation filtering, belonging to the field of point cloud extraction technology. Specifically, it includes: a point cloud preprocessing module for denoising and coordinate transformation of the original point cloud; an adaptive block division module for dividing the point cloud into uniform data blocks; a multi-scale feature extraction module for extracting elevation variation features, point cloud distribution density features, and surface texture features of each point under multiple neighborhood windows, forming a multi-scale visual cognitive feature set; a dynamic filtering module for normalizing and weighted fusion of features to generate comprehensive feature values, and for filtering candidate ground points using dynamic thresholds; a region growing and aggregation module for iteratively growing and aggregating candidate points based on feature similarity to form an optimized ground point set; and a boundary fusion output module for verifying feature consistency and adjusting coordinates at the data block splicing boundaries, outputting the final ground points. This invention achieves high-precision, adaptive extraction of LiDAR ground points in complex scenes.

Description

Technical Field

[0001] This invention relates to the field of point cloud extraction technology, and more specifically to a lidar ground point extraction system based on visual cognitive analog filtering. Background Technology

[0002] LiDAR, as an advanced 3D spatial information sensing technology, can quickly and accurately acquire discrete 3D point cloud data of a target scene by emitting laser beams and receiving reflected signals. This data plays an irreplaceable role in fields such as autonomous driving environmental perception, digital terrain model construction, urban 3D modeling, and natural resource surveys. Ground point extraction is a crucial preliminary step in LiDAR point cloud data processing. Its purpose is to accurately and completely separate the ground point set representing the real terrain surface from a mixed point cloud containing various terrain features such as ground, vegetation, buildings, and vehicles, providing reliable basic data for subsequent navigation, mapping, and analysis.

[0003] Currently, common LiDAR ground point extraction methods can be broadly categorized into two types. One type is based on mathematical models. These methods typically presuppose that the ground surface conforms to a specific geometric or mathematical model, such as progressive densification of triangular meshes (TINs) or planar / surface fitting, and then identify ground points by judging the degree of conformity between the point cloud and the presupposed model. The other type is based on filtering. These methods primarily design filtering rules or algorithms based on the spatial distribution characteristics of the point cloud (such as elevation, slope, and density), such as morphological opening operations and slope filtering, and distinguish between ground and non-ground points through iteration or threshold comparison. These methods can achieve certain results in relatively flat and simple scenes.

[0004] However, the aforementioned existing technologies still face significant challenges in dealing with complex real-world scenarios. Methods based on mathematical models heavily rely on prior terrain assumptions. In areas with significant terrain undulations or complex structures, such as steep slopes, gullies, vegetated areas, or densely built-up areas, the pre-set models often fail to accurately describe the true ground morphology, leading to missed ground point extraction or incorrect extraction of non-ground points. On the other hand, filtering methods based on fixed parameters or rules lack adaptability. Their key parameters, such as filter window size and slope threshold, are usually globally fixed values, unable to adapt to spatial variations in point cloud density, terrain undulation, and the complexity of ground features within the scene. In highly heterogeneous scenes, this can easily result in loss of detail or insufficient filtering. Therefore, there is an urgent need for a ground point extraction technology that can adapt to scene characteristics, does not rely on fixed prior models, and can robustly handle complex terrain. Summary of the Invention

[0005] The purpose of this invention is to provide a lidar ground point extraction system based on visual cognitive analog filtering, and to solve the following technical problems: Existing mathematical model-based methods rely on prior terrain assumptions, making it difficult to adapt to terrain changes such as steep slopes and vegetated areas, leading to misjudgments; while methods based on fixed parameters cannot dynamically adjust according to point cloud density and terrain undulation, and are prone to losing details or insufficient filtering in heterogeneous scenes.

[0006] The objective of this invention can be achieved through the following technical solutions: A lidar ground point extraction system based on visual cognitive analog filtering includes: The point cloud preprocessing module is used to input the raw point cloud data collected by the lidar, remove noise points and unify coordinate transformation of the raw point cloud data to obtain a preprocessed point cloud data set. The adaptive block segmentation module is used to divide the preprocessed point cloud data into multiple spatially continuous and non-overlapping point cloud data blocks according to the horizontal distribution range of the preprocessed point cloud data. The number of points contained in each point cloud data block is within a preset range. The multi-scale feature extraction module is used to extract the elevation change features, point cloud distribution density features, and surface texture features of each point in a neighborhood window of multiple different sizes for each point cloud data block, forming a multi-scale visual cognition feature set. The dynamic filtering module is used to normalize and weighted fuse the multi-scale visual cognitive features of each point in the point cloud data block to generate a comprehensive feature value. Points whose comprehensive feature values ​​fall within the dynamic threshold range are marked as candidate ground points, forming a set of candidate ground points. The region growth and aggregation module is used to select points whose feature vectors are located in the statistical center region from the candidate ground point set as the growth starting point, calculate the multi-scale feature similarity between the starting point and the neighboring points, and perform iterative growth and aggregation based on the similarity judgment results to form an optimized ground point set. The boundary fusion output module is used to spatially stitch together the optimized ground point set obtained after processing all point cloud data blocks, perform feature consistency verification and position adjustment on the points at the stitching boundary, and output the final ground point extraction result.

[0007] As a further aspect of the present invention: the process of dividing point cloud data into blocks in the adaptive block segmentation module is as follows: Obtain the two-dimensional spatial range of the preprocessed point cloud data projected onto the horizontal plane, calculate the minimum bounding rectangle of the two-dimensional spatial range, and construct a uniform two-dimensional grid with the lower left corner of the minimum bounding rectangle as the origin of the coordinate system. The two-dimensional grid completely covers the minimum bounding rectangle. Count the number of point clouds falling into each initial grid cell. For grid cells with a point cloud number exceeding a preset upper limit, divide them equally along their row and column directions to generate four sub-cells. For grid cells with a point cloud number below a preset lower limit, find the neighboring cell with the fewest point cloud numbers that is adjacent to the edge of the grid cell and merge it. The grid cells are iteratively divided and merged until the number of point clouds in all grid cells is within a preset range. Each final grid cell that meets the point cloud number range and its internal point cloud data are defined as an independent point cloud data block.

[0008] As a further aspect of the present invention: the construction process of multiple neighborhood windows of different sizes in the multi-scale feature extraction module is as follows: Calculate the average number of points per unit area within the point cloud data block, define the average number of points as the base density of the point cloud data block, calculate the reciprocal of the base density and take the square root, and use the calculation result as the reference window size. Set a fixed scaling factor greater than 1 and a positive integer representing the number of windows. Using the reference window size as the starting size, multiply the starting size by the fixed scaling factor to obtain the next window size. Repeat this process until the number of square windows with progressively increasing sizes is generated. The size of the last window does not exceed a preset ratio of the length of the diagonal of the minimum outer rectangle of the point cloud data block on the horizontal plane. Arrange the generated square windows in ascending order of size to form a multi-scale neighborhood window set for the point cloud data block.

[0009] As a further aspect of the present invention: the process of extracting elevation variation features, point cloud distribution density features, and surface texture features in the multi-scale feature extraction module is as follows: For a point to be processed within a point cloud data block, the center of each square window in the multi-scale neighborhood window set is aligned with the horizontal coordinate of the point to be processed. For each aligned square window, all point cloud data within its boundary are collected to form a neighborhood point set. The arithmetic mean of the elevation values ​​of all points in the neighborhood point set and the standard deviation of these elevation values ​​are calculated. The combination of the arithmetic mean and the standard deviation is recorded as the elevation change feature of the point to be processed at the current window scale. The total number of point clouds in the neighborhood point set is counted, and the total number is divided by the area of ​​the current square window. The result is recorded as the point cloud distribution density feature. The neighborhood point set is projected onto a horizontal plane, and a corresponding two-dimensional grayscale image is generated using the elevation value of the point as the grayscale value. The grayscale co-occurrence matrix of the two-dimensional grayscale image is calculated, and four statistical quantities, energy, contrast, correlation and entropy, are extracted from it. The numerical values ​​of these four statistical quantities are combined and recorded as surface texture features.

[0010] As a further aspect of the present invention: the process of generating comprehensive feature values ​​in the dynamic filtering module is as follows: For a multi-scale visual cognitive feature set of a point within a point cloud data block, traverse each feature value in the set, find the maximum and minimum values ​​of all points on the same feature dimension, and linearly scale the feature value of the current point using the maximum and minimum values ​​to map it to the range of zero to one. Calculate the standard deviation of the elevation values ​​of all points within the point cloud data block, and define the standard deviation as the terrain undulation of the data block. Assign a weight coefficient that is directly proportional to the terrain undulation to the elevation variation feature, assign a weight coefficient that is inversely proportional to the terrain undulation to the surface texture feature, and assign a fixed weight coefficient to the point cloud distribution density feature. Multiply each normalized feature value of the current point by its corresponding weight coefficient and sum them up. The sum is the comprehensive feature value of the point.

[0011] As a further aspect of the present invention: the generation process of the dynamic threshold interval in the dynamic filtering module is as follows: The comprehensive feature values ​​of all points within the point cloud data block are statistically analyzed, and the arithmetic mean and standard deviation of these comprehensive feature values ​​are calculated. The point cloud data block is uniformly divided into multiple rectangular sub-regions on a horizontal plane, and the point cloud density in each sub-region is calculated. The variance of the point cloud density in all sub-regions is calculated, and the variance is defined as the density spatial heterogeneity of the point cloud data block. A preset basic adjustment coefficient is obtained, and the basic adjustment coefficient is multiplied by the density spatial heterogeneity to obtain the dynamic adjustment coefficient. The product of the standard deviation and the dynamic adjustment coefficient is subtracted from the arithmetic mean of the comprehensive feature values, and the result is used as the lower limit of the dynamic threshold interval. The product of the standard deviation and the dynamic adjustment coefficient is added to the arithmetic mean, and the result is used as the upper limit of the dynamic threshold interval.

[0012] As a further aspect of the present invention: in the region growing and aggregation module, the process of calculating multi-scale feature similarity is as follows: Obtain the multi-scale visual cognitive feature set for each point in the candidate ground point set. For any two points, calculate the sum of the absolute values ​​of the differences between their corresponding feature values ​​as the feature difference degree. For each point in the set, calculate its feature difference degree with all other points in the set and sum them to obtain the overall difference degree of that point. Sort the candidate ground points in order of increasing overall difference degree. Select the first few points from the sorted sequence to form the initial growth starting point set. For each point in the initial growth starting point set, calculate its feature difference degree with all other points in the candidate ground point set. This feature difference degree is used as the inverse vector degree of the multi-scale feature similarity between the two points.

[0013] As a further aspect of the present invention: in the region growth and aggregation module, the process of iterative growth and aggregation based on the similarity judgment result to form an optimized set of ground points is as follows: Create a global set of processed point labels. Select an unlabeled point from the initial growth starting point set as the current growth seed point. Create an empty current growth subset and add the current growth seed point to it. Obtain the feature difference between the current growth seed point and all unlabeled candidate ground points. Add points with feature difference less than the preset growth similarity threshold to the current growth subset and record these new points as newly added growth seed points in this round. Mark all points in the current growth subset as processed. For each newly added seed point, the process of obtaining feature difference, determining the inclusion condition, and updating the subset and label is repeated until no new points are added. At this point, the current growth subset constitutes an optimized ground point subset. The process of selecting seed points from the initial growth starting point set, creating and growing subsets is repeated until all starting points are marked as processed. Finally, the union of all subsets constitutes the optimized ground point set.

[0014] As a further aspect of the present invention: the process of outputting the final ground point extraction result in the boundary fusion output module is as follows: Merge the optimized ground point subsets of all point cloud data blocks to form an initial global point cloud. Set a boundary buffer distance and identify all points whose horizontal projected coordinates are less than the boundary buffer distance of any point cloud data block boundary line. Mark them as the boundary point set to be verified. For each point in the boundary point set to be verified, find its neighboring points from different point cloud data blocks and calculate the Euclidean distance between the feature set of the point and the average feature set of the neighboring points. The Euclidean distance is compared with a preset consistency threshold. If the threshold is exceeded, the original feature set of the point is replaced with the average feature value of the neighboring points. For all boundary points whose feature sets have been replaced, their three-dimensional coordinates are recalculated by distance-weighted interpolation based on their new feature values ​​and the coordinates of the neighboring points. After all boundary points are verified and their coordinates are adjusted, the final ground point extraction result is output.

[0015] The beneficial effects of this invention are: This invention effectively overcomes the dependence of existing technologies on fixed mathematical models and global parameters by simulating the multi-scale perception process of visual cognition. Specifically, it employs an adaptive multi-scale neighborhood window for feature extraction, which can simultaneously capture local details and global trends, thus maintaining the integrity of terrain features even in steep slopes, gullies, and areas with mixed vegetation. By introducing a dynamic weight allocation and threshold generation mechanism based on the spatial heterogeneity of terrain undulation and density in data blocks, it achieves complete adaptation of filtering standards, significantly improving the robustness and generalization ability of the method under different density and undulation scenarios. Furthermore, a region growing algorithm based on multi-scale feature similarity can accurately aggregate points based on their own features, effectively distinguishing non-ground points with similar appearances. Finally, by verifying the consistency of boundary features and adjusting coordinate interpolation, it eliminates the seams introduced by block processing, ensuring the overall smoothness and spatial continuity of the output ground point set. The entire scheme forms a closed loop from feature perception and decision-making to result optimization, achieving high-precision and high-completeness terrain extraction in complex and heterogeneous scenes without the need for a pre-set ground model. Attached Figure Description

[0016] The invention will now be further described with reference to the accompanying drawings.

[0017] Figure 1 This is a schematic diagram of the process of a lidar ground point extraction system based on visual cognitive analog filtering according to the present invention. Detailed Implementation

[0018] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0019] Please see Figure 1 As shown, this invention is a lidar ground point extraction system based on visual cognitive analog filtering, comprising: The point cloud preprocessing module is used to input the raw point cloud data collected by the lidar, remove noise points and unify coordinate transformation of the raw point cloud data to obtain a preprocessed point cloud data set. The adaptive block segmentation module is used to divide the preprocessed point cloud data into multiple spatially continuous and non-overlapping point cloud data blocks according to the horizontal distribution range of the preprocessed point cloud data. The number of points contained in each point cloud data block is within a preset range. The multi-scale feature extraction module is used to extract the elevation change features, point cloud distribution density features, and surface texture features of each point in a neighborhood window of multiple different sizes for each point cloud data block, forming a multi-scale visual cognition feature set. The dynamic filtering module is used to normalize and weighted fuse the multi-scale visual cognitive features of each point in the point cloud data block to generate a comprehensive feature value. Points whose comprehensive feature values ​​fall within the dynamic threshold range are marked as candidate ground points, forming a set of candidate ground points. The region growth and aggregation module is used to select points whose feature vectors are located in the statistical center region from the candidate ground point set as the growth starting point, calculate the multi-scale feature similarity between the starting point and the neighboring points, and perform iterative growth and aggregation based on the similarity judgment results to form an optimized ground point set. The boundary fusion output module is used to spatially stitch together the optimized ground point set obtained after processing all point cloud data blocks, perform feature consistency verification and position adjustment on the points at the stitching boundary, and output the final ground point extraction result.

[0020] In a preferred embodiment of the present invention, the specific working process of the adaptive block segmentation module is as follows: Obtain the 2D projection range of the preprocessed point cloud data onto the horizontal plane. This range is determined by calculating the minimum and maximum X-coordinates (Xmin and Xmax) and the minimum and maximum Y-coordinates (Ymin and Ymax) of all points. Therefore, define a rectangular region that completely encloses the horizontal projection of all points; this rectangle is called the minimum bounding rectangle, with its lower left vertex coordinates (Xmin, Ymin) and upper right vertex coordinates (Xmax, Ymax).

[0021] Next, using the bottom left vertex of the smallest bounding rectangle as the origin of the coordinate system, an initial two-dimensional uniform grid covering the entire rectangle is constructed. The size of the initial grid, i.e., the side length of each cell, can be set based on the total number of points in the point cloud and empirical values ​​for the application scenario. For example, for a medium-density area of ​​approximately one square kilometer, the initial grid side length can be set to 50 meters. This grid consists of square cells arranged in a regular pattern of multiple rows and columns, ensuring that each point in the original point cloud can be assigned to a unique initial grid cell based on its X and Y coordinates.

[0022] The system then iterates through all initial grid cells, counting the number of point clouds falling into each cell. This step aims to evaluate the uniformity of the point cloud distribution on the horizontal plane. To ensure the efficiency and stability of subsequent processing, the number of point clouds in each data processing unit needs to be controlled within a reasonable range, such as between 2000 and 3000 points. This range is the preset interval, with its lower limit called the preset lower limit and its upper limit called the preset upper limit.

[0023] Based on the above statistics, the system dynamically optimizes the grid structure, with core operations including the equal division and merging of grid cells.

[0024] For cells with a point cloud count exceeding the preset limit, it indicates that the point cloud in that area is too dense and needs to be refined. The equal division operation is performed simultaneously along both the row direction (X-axis) and column direction (Y-axis) of the cell. Specifically, find the center point of the cell, and draw two straight lines parallel to the X-axis and Y-axis through that center point, dividing the original cell into four equal-area sub-cells. Then, recount the point cloud count of each of these four sub-cells.

[0025] For cells with a point cloud count below a preset lower limit, it indicates that the point cloud in that area is too sparse and needs to be merged with adjacent areas to form a processing unit with a suitable data volume. The merging operation first locates all edge-adjacent cells of the current cell, i.e., cells that share a complete edge with it. Then, among these adjacent cells, the one with the fewest point cloud counts is selected. Finally, the current cell is merged with this selected adjacent cell to form a larger polygonal cell, which may be rectangular in shape. After merging, the point cloud count of the new cell is the sum of the point cloud counts of the original two cells.

[0026] The above division and merging operations are not performed once and then stop; they are iterative. After a complete loop through all cells, the system checks the point cloud count in each existing cell again. As long as there are still cells whose point cloud count is outside the preset range, a new round of traversal and adjustment operations is initiated. This is a process of repeated iterations until the conditions are met. The criterion for terminating the iteration is: among all currently existing grid cells, no cell has a point cloud count below the preset lower limit, and no cell has a point cloud count above the preset upper limit.

[0027] Once the iteration process is complete, an adaptively optimized mesh partitioning scheme is obtained. At this point, each final mesh cell satisfies the condition that the number of point clouds is within a preset range. The three-dimensional space covered by such a final mesh cell, along with all point cloud data falling within this space, is defined as an independent basic unit for subsequent processing, called a point cloud data block.

[0028] To facilitate management and tracking, the system assigns a unique spatial identifier to each point cloud data block upon generation. The encoding rules for this identifier are typically strongly correlated with its spatial location. One feasible approach is to use the coordinates of the bottom-left vertex of the final grid cell corresponding to the data block for encoding. For example, for a data block with bottom-left vertex coordinates of 100.0, 200.0, its identifier could be generated as Block100.0 200.0. This identifier will play a crucial role in the subsequent data stitching and result integration stages, accurately tracing the original processing unit for each point.

[0029] In another preferred embodiment of the present invention, the specific working process of the multi-scale feature extraction module is as follows: Phase 1: Construction of Multi-Scale Neighborhood Windows First, a base window size needs to be determined that matches the characteristics of the current point cloud data block. This size depends on the point cloud density of the data block. The average number of points per unit area within the point cloud data block is calculated; this value is defined as the "base density" of the data block. For example, assuming a point cloud data block contains 10,000 points and its projected area on a horizontal plane is 2,500 square meters, its base density is 4 points per square meter. The base density reflects the average sparsity of the point cloud in that area.

[0030] Using the base density directly as a reference for window size is not intuitive because density is measured in "points per area," while window size is a unit of length. Therefore, a mathematical transformation is needed to establish the connection. A reasonable and effective transformation is to calculate the reciprocal of the base density and then take its square root. Physically, this calculation approximates the expected neighborhood radius that can contain a constant number of point clouds at the average density. Continuing with the example above, with a base density of 4 points per square meter, its reciprocal is 0.25 square meters per point, and taking its square root yields 0.5 meters. This 0.5 meters is the calculated "baseline window size." It means that at the current average density, a square area with sides of approximately 0.5 meters is expected to contain about one point. This provides a reasonable starting size for the finest-scale analysis.

[0031] After obtaining the baseline window size, a series of windows with progressively increasing sizes can be constructed. This process requires two preset parameters: a "fixed scaling factor" greater than 1 and a "positive integer" representing the total number of windows. The fixed scaling factor determines the size increase between adjacent scale windows and is usually set according to the desired scale range, such as 1.5, 1.8, or 2.0. The number of windows determines the number of scales, such as 3, 4, or 5, to balance computational complexity and feature richness.

[0032] When constructing the sequence, the baseline window size is used as the first window, which is also the side length of the smallest window. The size of the first window is multiplied by a fixed scaling factor to obtain the size of the second window. The size of the second window is then multiplied by the same fixed scaling factor to obtain the size of the third window. This multiplication process is repeated until a preset number of windows are generated. For example, if the baseline size is 0.5 meters, the fixed scaling factor is 2.0, and the number of windows is 4, then the generated window size sequence is: 0.5 meters, 1.0 meters, 2.0 meters, and 4.0 meters. This forms a set of square windows with sizes increasing geometrically.

[0033] To prevent the maximum window size from exceeding the effective analysis range of the current data block, an upper limit needs to be set. Typically, the side length of the last window is limited to a certain proportion of the diagonal length of the minimum horizontal bounding rectangle of the point cloud data block, such as one-half or one-third. If the calculated maximum window size exceeds this limit, the windows exceeding the limit are discarded, or the fixed proportion coefficient and the number of windows are adjusted and regenerated. Finally, all the generated square windows are arranged in ascending order of their side lengths, forming a "multi-scale neighborhood window set" for the point cloud data block. This set allows for subsequent observation and analysis of the same point at different spatial scales.

[0034] Phase Two: Extraction of Multi-Scale Visual Cognitive Features For each point to be processed within the point cloud data block, the multi-scale neighborhood window set constructed above will be used to sequentially calculate its features in three different dimensions: elevation variation features, point cloud distribution density features, and surface texture features.

[0035] Step 1: Obtaining the Neighborhood Point Set. For the point to be processed, obtain its X and Y coordinates on the horizontal plane. Then, sequentially select a window from the multi-scale neighborhood window set. Using the X and Y coordinates of the point to be processed as the geometric center of the current square window, determine the precise boundary of the window on the horizontal plane. Subsequently, search and collect all points whose horizontal coordinates fall within the boundary of this square window. These collected points constitute the "neighborhood point set" of the current point to be processed at the current specific window scale. It is important to emphasize that the neighborhood determination here is based solely on the horizontal coordinates; the elevation Z coordinate is not involved in the delineation of the neighborhood range.

[0036] Step Two: Calculation of Elevation Variation Characteristics. After obtaining the neighborhood point set, first calculate the arithmetic mean of the elevation Z values ​​of all points. This mean represents the average elevation level within the neighborhood. Next, calculate the standard deviation of these elevation values ​​relative to this mean. The standard deviation quantifies the elevation fluctuation or dispersion of points within the neighborhood. For example, on flat ground, the elevation standard deviation will be small; while on a thicket or small mound, the elevation standard deviation will be significantly larger. The calculated arithmetic mean and standard deviation are recorded as a numerical pair, collectively as the "elevation variation characteristic" of the point to be processed at this window scale. This characteristic includes both absolute positional information and relative change information of the elevation.

[0037] Step 3: Calculation of Point Cloud Distribution Density Features. This feature aims to describe the local aggregation of point clouds on a horizontal plane. The calculation process is relatively straightforward: count the total number of point clouds contained in the current neighborhood point set. Then, divide this total number by the area of ​​the currently used square window (i.e., the square of the window's side length). The result is a density value, expressed in points per square meter. This value is recorded as the "point cloud distribution density feature" of the point to be processed at this window scale. For example, on a solid, flat surface, laser reflection is stable, and the point cloud density may be high and uniform; while on sparse grass, the point cloud density may be lower and more variable.

[0038] Step 4: Surface Texture Feature Calculation. This feature aims to transform the undulation patterns of the local surface of a 3D point cloud into quantifiable image texture indicators, simulating the visual perception of surface roughness and regularity. First, the current neighborhood point set is orthogonally projected onto a horizontal plane. After projection, each point retains only its X,Y coordinates, but its elevation Z value is not discarded; instead, it is used as a key attribute. Next, a discrete pixel grid is established within the 2D region covered by the current square window. The grid resolution needs to be set, for example, one pixel per 0.1 meters. For each pixel, all points that fall within the projection range are found. If multiple points exist, the average elevation Z value of these points is taken as the grayscale value of the pixel; if no points fall within the range, interpolation or a background value needs to be assigned. Through this process, a 2D grayscale image is generated, where the brightness (grayscale) of a pixel directly corresponds to the elevation of the point cloud at that location.

[0039] Next, texture analysis is performed on the generated grayscale image using the classic tool, the "Gray-Level Co-occurrence Matrix" (GLCM). The GLCM describes the frequency of pixel pairs with a certain grayscale relationship (such as specific distance and orientation) in an image. First, the parameters for calculating this matrix need to be set: one is the "grayscale level," which determines how many levels the continuous grayscale values ​​are quantized to, for example, 16 levels; the second is the "distance," which determines the distance between pixel pairs to consider, for example, a distance of 1 pixel; and the third is the "orientation," typically considering four main directions: 0 degrees, 45 degrees, 90 degrees, and 135 degrees.

[0040] During calculation, all pixel pairs in the image are traversed according to a set distance and direction. The number of occurrences of pixel pairs with different grayscale combinations is counted, forming a square matrix. To obtain a rotation-invariant texture description, the matrices calculated in four different directions are usually averaged to obtain a comprehensive grayscale co-occurrence matrix.

[0041] Finally, extract the four most commonly used statistics from this composite matrix: Energy: Reflects the uniformity of gray-level distribution and the coarseness of texture in an image. A higher energy value indicates a highly concentrated gray-level distribution. Contrast: Reflects the sharpness of an image and the depth of texture grooves. Greater local gray-level variations in an image result in a higher contrast value. Correlation: Measures the degree of linear correlation between gray levels in an image along rows or columns, reflecting the directionality of texture. Entropy: Represents the amount of information in an image or the randomness of texture. More complex and random textures result in a higher entropy value.

[0042] The calculated energy, contrast, correlation, and entropy values ​​are combined into an array and recorded as the "surface texture features" of the point to be processed at this window scale.

[0043] Step 5: Feature Set Integration. For the same point to be processed, repeat steps one through four, traversing all windows in the multi-scale neighborhood window set. Ultimately, a set of feature values ​​for the point at each window scale will be obtained, including an elevation variation feature pair, a point cloud distribution density value, and a surface texture feature array containing four components. All these feature values ​​are concatenated into a one-dimensional long vector in ascending order of window scale. This long vector is the point's "multi-scale visual cognitive feature set." It comprehensively encodes the geometric, statistical, and textural attributes of the point and its surrounding local area in three-dimensional space from multiple scales and dimensions, providing a rich and robust data foundation for subsequent intelligent filtering and classification.

[0044] In another preferred embodiment of the present invention, the specific working process of the dynamic filtering module is as follows: Phase 1: Generation of Comprehensive Feature Values First, feature normalization is performed. For a point cloud data block, assume it contains N points, each with a multi-scale visual cognitive feature set of dimension D. Since these features originate from different physical quantities, such as the standard deviation of elevation, point cloud density, and texture energy, their numerical ranges and units differ. To consider each feature dimension equally in subsequent fusion, each feature dimension needs to be normalized independently. Specifically, for the d-th feature dimension, iterate through all N points in the data block and find the maximum value Maxd and the minimum value Mind among the feature values ​​of that dimension. These two values ​​define the numerical distribution range of the current data block in that feature dimension. For any point Pi in the data block, its original feature value Vi,d in the d-th dimension is mapped to a closed interval between zero and one through the following linear transformation: calculate the normalized value V'i,d = (Vi,d - Mind) / (Maxd - Mind). When Maxd equals Mind, i.e., all values ​​in that dimension are the same, the normalization value is defined as 0.5. By performing this operation on each of the D feature dimensions, the entire feature set of point Pi is transformed into a normalized feature vector where all elements are between 0 and 1. This step eliminates the influence of dimensions, making different features comparable.

[0045] Secondly, scene-adaptive weight allocation is performed. Simple weighted summation is not the essence of this method; the determination of weights needs to reflect the terrain characteristics of the current scene. Therefore, a key indicator describing the scene's topography, "terrain relief," is introduced. The calculation of terrain relief is straightforward: extract the elevation values ​​of all points within the point cloud data block and calculate the standard deviation of these elevation values. For example, the standard deviation of elevation for a flat parking lot might be only 0.05 meters; while for a grassland containing slopes and hills, the standard deviation might reach 0.5 meters or higher. This standard deviation, i.e., terrain relief, quantifies the overall degree of vertical fluctuation of the data block.

[0046] The weighting strategy is closely related to the terrain undulation, and its core logic is to simulate the degree to which human vision relies on various cues under different terrain conditions. The specific rules are as follows: Weighting of Elevation Variation Features: This weight is set to be proportional to the terrain undulation. The principle is that in areas with dramatic terrain undulation, the elevation variation of the ground itself is the most critical and stable feature distinguishing it from upright objects such as buildings. The elevation variation features of a point on a steep slope will be significantly different from those of non-ground vegetation clusters. Therefore, when the terrain undulation is high, the elevation variation feature should be given a higher weight to enhance its ability to identify ground points. A simple implementation is to set its weighting coefficient Wz = a × TerrainRoughness + b, where a is a positive constant, b is a base value, and TerrainRoughness is the calculated terrain undulation.

[0047] Weighting of surface texture features: This weight is set to be inversely proportional to the terrain undulation. This is because, in very flat areas, the distinguishability of elevation variation features and point cloud density features decreases. For example, a flat asphalt road surface and the roof of a truck may have similar low elevation variation and similar point cloud density. In this case, the surface texture pattern becomes the key discrimination criterion. The texture of a smooth road surface is uniform and regular, while the texture pattern of a truck roof may vary due to structural components. Therefore, when the terrain undulation is low, texture features should be given a higher weight. Its weighting coefficient Wt = c / (TerrainRoughness + e) ​​+ d, where c, d, and e are constants, and e is used to prevent division by zero.

[0048] Weighting of point cloud distribution density features: This feature is considered to have relatively stable distinguishing value under different terrains, and therefore is assigned a fixed weighting coefficient Wρ. This fixed value can be set empirically.

[0049] Finally, weighted fusion is performed. For point Pi, each value in its normalized feature vector is multiplied by the weight coefficient corresponding to the category of that feature dimension, calculated in the previous steps. Then, all weighted feature values ​​are summed. This final sum is the "comprehensive feature value" of point Pi. This scalar value encapsulates the point's location information in a multi-scale, multi-dimensional feature space, and its calculation process has been adjusted for importance based on the current local terrain.

[0050] Phase Two: Generation of Dynamic Threshold Ranges The purpose of generating a dynamic threshold range is to set an adaptive and reasonable judgment threshold for the comprehensive feature value to filter candidate ground points. The core idea is that the threshold should be determined based on the overall performance of all points within the current data block, and its tolerance should be able to adapt to the non-uniformity of the point cloud distribution.

[0051] First, calculate the global statistics. For a data block containing N points, the composite eigenvalue Fi for each point has already been calculated. Calculate the arithmetic mean μF of these N composite eigenvalues, and their standard deviation σF. The mean μF represents the “average eigenvalue level” of all points within the data block, while the standard deviation σF reflects how much these points are dispersed around the average value.

[0052] Secondly, the uniformity of the point cloud distribution is assessed, i.e., the "density spatial heterogeneity" is calculated. The degree of spatial variation in point cloud density directly affects the local statistical validity of features. To quantify this variation, the entire point cloud data block is uniformly divided into M rectangular sub-regions on a horizontal plane. For example, a 100-meter square data block can be divided into 10 rows and 10 columns, totaling 100 sub-regions. Then, for each sub-region j, the point cloud density ρj within it is calculated. Next, the variance Varρ of these M sub-region density values ​​ρ1, ρ2, ..., ρM is calculated. This variance Varρ is defined as the "density spatial heterogeneity" of the point cloud data block. If the point cloud distribution is very uniform, and the density of each sub-region is similar, the variance is small, and the heterogeneity is low; if the point cloud is dense in some areas and sparse in others, the variance is large, and the heterogeneity is high. High heterogeneity means that multiple ground features or areas with large density differences may coexist within the scene, requiring a more tolerant threshold to handle this uncertainty.

[0053] Next, the dynamic adjustment coefficient is calculated. A preset "base adjustment coefficient" Kbase is obtained, which represents the threshold relaxation baseline under an ideal uniform distribution. Then, the base adjustment coefficient Kbase is multiplied by the calculated density spatial heterogeneity Varρ to obtain the final "dynamic adjustment coefficient" Kdynamic, i.e., Kdynamic = Kbase × Varρ. When the heterogeneity Varρ is high, Kdynamic increases; when the heterogeneity is low, Kdynamic approaches Kbase.

[0054] Finally, the dynamic threshold interval is determined. Using the previously calculated average μF, standard deviation σF, and dynamic adjustment coefficient Kdynamic of the global integrated eigenvalues, the threshold interval is determined using the following formula: The lower threshold Tlow = μF - σF × Kdynamic, and the upper threshold Thigh = μF + σF × Kdynamic.

[0055] The resulting interval [Tlow, Thigh] is a dynamic range centered on the global average feature level, with a width proportional to the feature dispersion and density heterogeneity. The more uneven the density distribution and the more dispersed the feature values, the wider the threshold interval should be to accommodate greater internal variability. After generating the threshold interval, all points within the data block are traversed, and points whose comprehensive feature value Fi falls within this interval are marked as "candidate ground points," thus forming a set of candidate ground points. This process achieves adaptive filtering that is entirely data-driven and does not require a preset fixed threshold.

[0056] In another preferred embodiment of the present invention, the specific working process of the region growth and aggregation module is as follows: Phase 1: Quantification of multi-scale feature similarity and selection of initial seed points The quality of region growing is highly dependent on the accurate measurement of similarity between points and the appropriate selection of the starting point for growth. This method uses numerical differences based on multi-scale visual cognitive feature sets to define similarity and selects the initial growth core based on the global representativeness of the points.

[0057] First, feature dissimilarity is calculated. Assume that after dynamic filtering, a point cloud data block yields a candidate ground point set containing M points. Each point Pi possesses a multi-scale visual cognitive feature set of dimension D, which is a comprehensive vector containing elevation, density, and texture information at multiple scales. To measure the similarity between any two points Pi and Pj, a computable "feature dissimilarity" metric needs to be defined. This metric is calculated directly and effectively: for each dimension in the feature sets of the two points, the absolute value of the difference between the corresponding feature values ​​is calculated. Then, the absolute differences across all D dimensions are summed. This can be expressed by the formula: Feature Dissimilarity Diiffij = Σ{k=1 to D}|Vi,k-Vj,k|, where Vi,k represents the feature value of point Pi in the k-th dimension. The smaller this Diiffij value, the closer the two points are across all feature dimensions, i.e., the more similar they are; conversely, the larger the Diiffij value, the greater the difference between the two points. This method comprehensively considers the differences in point cloud attributes in terms of geometry, statistics, and texture.

[0058] Secondly, the overall dissimilarity of each point is calculated and ranked. To evaluate the overall consistency or representativeness of a single point with the entire candidate ground point set, it is necessary to calculate the "overall dissimilarity" of that point relative to all other points in the set. For point Pi, its overall dissimilarity TotalDiffi is obtained by calculating its feature dissimilarity Diiffij with every other point Pj in the set, and summing all these Diiffij values, i.e., TotalDiffi = Σ{j=1,j≠i to M}Diffij. This TotalDiffi value has a clear physical meaning: it quantifies the degree to which the features of point Pi deviate from the average features of all other points in the set. A point with a small TotalDiffi value means that its features are very similar to the vast majority of points in the set, and it is likely located in the core region of the feature distribution, making it a typical and representative ground point. Conversely, a point with a large TotalDiffi value may be an outlier with singular features, or a boundary point located in the transition region between the ground and non-ground areas.

[0059] Next, all M candidate ground points are sorted from smallest to largest based on their overall dissimilarity. After sorting, the point at the very beginning of the sequence, i.e., the point with the smallest overall dissimilarity, is considered the most "ordinary" and representative point. These points are ideal cores for initiating region growth because they are located in dense areas of the feature space, and growing around them makes it easier to aggregate a large number of similar points.

[0060] Finally, an initial set of growth starting points is selected. From the sorted sequence described above, the first N points are chosen to form the "initial growth starting point set." The number N can be set as a fixed percentage of the total number of points M, such as 5% to 10%, or it can be set as an absolute value, such as 50 points. The key is that these selected points are the core group of points with the highest feature consistency; they are scattered throughout the data block, providing multiple reliable starting points for subsequent parallel region growth. After selection, it is necessary to pre-calculate the feature difference (Diff) between each point in the initial growth starting point set and all other points in the entire candidate ground point set, and store these calculation results for quick retrieval in the subsequent growth process. This data constitutes the direct basis for judging the similarity between two points.

[0061] Phase 2: Iterative growth and aggregation based on feature similarity This stage focuses on the initial growth starting point, and gradually expands to form a connected and optimized ground region by iteratively searching and merging points with similar features.

[0062] First, initialize the data structure. Create a global "processed point tag set", typically a boolean array or hash table of the same length as the candidate ground point set, to record whether each candidate ground point has been assigned to a certain grown subset. All points are initially labeled "unprocessed".

[0063] The growth process proceeds cyclically, processing each point in the initial set of growth starting points one by one. From this set, a point currently marked as "unprocessed" is selected sequentially and designated as the "current growth seed point." A new, empty set called the "current growth subset" is created for this seed point. The current growth seed point is added to this subset, and its status in the global label set is updated to "processed."

[0064] Next, a single growth iteration is performed with the current seed point as the core. The feature difference (Diff) between the current growth seed point and all candidate ground points marked "unprocessed" is queried (either pre-calculated or calculated in real-time). A "growth similarity threshold" (Thgrowth) is set. This threshold is a key parameter determining whether two points are similar enough to be merged; its value needs to be determined experimentally, for example, it can be set as a proportion of the average feature difference between all point pairs. All unprocessed points are traversed, and it is determined whether their Diff value with the current seed point is less than Thgrowth. If the condition is met, the point is considered sufficiently similar to the current seed point, added to the current growth subset, and simultaneously recorded as a "newly added growth seed point in this round." The global label of all newly added points is also updated to "processed." After completing one traversal of all unprocessed points, the current growth subset is initially expanded.

[0065] Then, the iterative expansion phase begins. The above process only performed a single growth step using an initial seed point. To fully aggregate all similar points, it is necessary to continue growth with newly added points as the new starting point. The list of "newly added seed points in this round" is used as the seed point set for the next round of growth. For each newly added seed point in this list, the aforementioned steps are repeated: calculate the feature difference between the newly added seed point and all remaining "unprocessed" points, add points with a difference less than Thgrowth to the current growth subset, and record these new points as newly added seed points for the next round. This process is continuously repeated, with each round expanding outwards using the newly added points from the previous round as seeds. The iteration terminates when no new points are added to the current growth subset in a complete traversal. This means that, starting from the current initial seed point, all points reachable within the feature similarity threshold have been aggregated. At this point, the current growth subset is closed, becoming an "optimized ground point subset" with highly consistent internal features and potentially spatial connectivity.

[0066] After processing an initial growth starting point, the program returns to the loop start point and searches for the next point still marked as "unprocessed" in the initial growth starting point set. Using this point as the new seed, the entire process of creating new subsets and iterative growth is repeated. In this way, each initial growth starting point can potentially grow an optimized subset of ground points. Because the growth process marks processed points, different subsets will not contain the same points.

[0067] Finally, when all points in the initial growth starting point set have been processed (i.e., have become seeds of a subset or have been incorporated into a subset), the entire growth and aggregation process ends. At this point, K optimized ground point subsets may have been generated, where K is less than or equal to the number of initial growth starting points N. The union of these K subsets, i.e., all successfully aggregated points, constitutes the "final optimized ground point set" of the point cloud data block. Candidate ground points that were not included in any subset during the entire process, even after passing dynamic filtering, are considered noise or misclassified points because their features differ too much from all the core growth points, and are thus removed at this stage. Through this iterative growth based on feature similarity, the method can effectively filter out discrete and inconsistent points while preserving the continuity of the main ground structure, thereby significantly improving the purity and regional consistency of the ground point extraction results.

[0068] In another preferred embodiment of the present invention, the specific working process of the boundary fusion output module is as follows: First, the optimized ground point subsets obtained after region growing and aggregating all point cloud data blocks are merged. The merging process is based on the spatial identifier code of the point cloud data block to which each point belongs and its original 3D coordinates, bringing all points together into a unified data structure to form the initial global point cloud. At this point, although the points at the boundaries of different data blocks are spatially close, their multi-scale visual cognitive feature sets may have slight inconsistencies because they are processed and judged in different local contexts.

[0069] To locate these areas requiring special attention, a boundary buffer distance needs to be defined. This distance defines a strip-shaped region on both sides of the data block boundary line. For example, if the common boundary of two adjacent point cloud data blocks is a straight line, then a buffer zone extending 2 meters on each side of this line constitutes the buffer zone. The program iterates through the boundary lines of all point cloud data blocks. For each point in the global point cloud, it calculates the vertical distance from its horizontal projected coordinates to the nearest data block boundary line. If this distance is less than the preset boundary buffer distance, such as 2 meters, then the point is marked as a boundary point to be verified, and all marked points constitute the set of boundary points to be verified.

[0070] Next, feature consistency verification is performed on each boundary point to be verified. Taking the current point P to be verified as an example, a search is performed within its 3D spatial neighborhood to find the K nearest neighbor points, where K is usually set to an integer between 10 and 20. The key to this search is that at least some of the found neighbor points must come from a different point cloud data block than point P. This can be achieved by checking the original data block identifier to which each neighbor point belongs. Only neighbor points from different data blocks can carry feature information that provides an independent reference for verification. The average of the multi-scale visual cognitive feature sets of these neighbor points from different data blocks is calculated to obtain a reference feature vector Vref. The reference feature vector Vref represents the feature consensus provided by other surrounding data blocks at the spatial location of point P. Subsequently, the Euclidean distance D between point P's own feature vector Vself and the reference feature vector Vref is calculated. The Euclidean distance D directly quantifies the degree of deviation between the features of point P and the feature consensus of the local region.

[0071] The calculated Euclidean distance D is compared with a preset consistency threshold Tcons. The consistency threshold Tcons is an empirical value determined based on a large amount of experimental data, which defines the upper limit of acceptable feature differences. If the distance D is greater than the threshold Tcons, it is determined that the features of point P are significantly inconsistent, and the original feature set Vself of point P needs to be replaced with a reference feature vector Vref to achieve feature smoothing and alignment.

[0072] For points whose feature sets are replaced, their spatial coordinates may also need fine-tuning to improve geometric continuity. Coordinate adjustment uses a distance-weighted interpolation algorithm. Taking point P as an example, the 3D coordinates of its neighboring points from different data blocks are used. Let the coordinates of the j-th neighboring point be Xj, Yj, Zj, and its Euclidean distance from the original coordinates of point P be dj. First, the weight wj of each neighboring point is calculated based on these distances. The weight is inversely proportional to the distance, and is usually calculated using wj = 1 / (dj + ε), where ε is a very small constant to prevent division by zero. Then, the new horizontal coordinates of point P are calculated using a weighted average: X'p = Σ(wj × Xj) / Σwj, Y'p = Σ(wj × Yj) / Σwj. The adjustment of the elevation coordinate Z'p requires incorporating its new feature vector. Since the feature vector contains elevation change information, after feature replacement, the corrected elevation value can be obtained by weighted interpolation based on the elevation statistics implied in the new feature vector and the elevations of neighboring points. Through this process, the position of the boundary point is moderately pulled toward the average position of points from other data blocks in its neighborhood, thereby reducing the gaps between blocks.

[0073] Finally, after traversing and processing all points in the set of boundary points to be verified, the features and coordinates of all points in the global point cloud are finally integrated. All the boundary points that have passed verification or have been adjusted, together with the interior points, constitute a final ground point extraction result dataset that is more continuous and consistent in terms of features and geometric space, and this dataset is used as the output of the entire method.

[0074] The foregoing has provided a detailed description of one embodiment of the present invention, but this description is merely a preferred embodiment and should not be construed as limiting the scope of the invention. All equivalent variations and modifications made within the scope of the claims of this invention should still fall within the patent coverage of this invention.

Claims

1. A laser radar ground point extraction system based on visual cognitive simulation filtering, characterized in that, include: The point cloud preprocessing module is used to input the raw point cloud data collected by the lidar, remove noise points and unify coordinate transformation of the raw point cloud data to obtain a preprocessed point cloud data set. The adaptive block segmentation module is used to divide the preprocessed point cloud data into multiple spatially continuous and non-overlapping point cloud data blocks according to the horizontal distribution range of the preprocessed point cloud data. The number of points contained in each point cloud data block is within a preset range. The multi-scale feature extraction module is used to extract the elevation change features, point cloud distribution density features, and surface texture features of each point in a neighborhood window of multiple different sizes for each point cloud data block, forming a multi-scale visual cognition feature set. The dynamic filtering module is used to normalize and weighted fuse the multi-scale visual cognitive features of each point in the point cloud data block to generate a comprehensive feature value. Points whose comprehensive feature values ​​fall within the dynamic threshold range are marked as candidate ground points, forming a set of candidate ground points. The region growth and aggregation module is used to select points whose feature vectors are located in the statistical center region from the candidate ground point set as the growth starting point, calculate the multi-scale feature similarity between the starting point and the neighboring points, and perform iterative growth and aggregation based on the similarity judgment results to form an optimized ground point set. The boundary fusion output module is used to spatially stitch together the optimized ground point set obtained after processing all point cloud data blocks, perform feature consistency verification and position adjustment on the points at the stitching boundary, and output the final ground point extraction result.

2. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 1, characterized in that, In the adaptive block segmentation module, the process of dividing the point cloud data into blocks is as follows: Obtain the two-dimensional spatial range of the preprocessed point cloud data projected onto the horizontal plane, calculate the minimum bounding rectangle of the two-dimensional spatial range, and construct a uniform two-dimensional grid with the lower left corner of the minimum bounding rectangle as the origin of the coordinate system. The two-dimensional grid completely covers the minimum bounding rectangle. Count the number of point clouds falling into each initial grid cell. For grid cells with a point cloud number exceeding a preset upper limit, divide them equally along their row and column directions to generate four sub-cells. For grid cells with a point cloud number below a preset lower limit, find the neighboring cell with the fewest point cloud numbers that is adjacent to the edge of the grid cell and merge it. The grid cells are iteratively divided and merged until the number of point clouds in all grid cells is within a preset range. Each final grid cell that meets the point cloud number range and its internal point cloud data are defined as an independent point cloud data block.

3. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 1, characterized in that, In the multi-scale feature extraction module, the construction process of multiple neighborhood windows of different sizes is as follows: Calculate the average number of points per unit area within the point cloud data block, define the average number of points as the base density of the point cloud data block, calculate the reciprocal of the base density and take the square root, and use the calculation result as the reference window size. Set a fixed scaling factor greater than 1 and a positive integer representing the number of windows. Using the reference window size as the starting size, multiply the starting size by the fixed scaling factor to obtain the next window size. Repeat this process until the number of square windows with progressively increasing sizes is generated. The size of the last window does not exceed a preset ratio of the length of the diagonal of the minimum outer rectangle of the point cloud data block on the horizontal plane. Arrange the generated square windows in ascending order of size to form a multi-scale neighborhood window set for the point cloud data block.

4. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 3, characterized in that, In the multi-scale feature extraction module, the process of extracting elevation change features, point cloud distribution density features, and surface texture features is as follows: For a point to be processed within a point cloud data block, the center of each square window in the multi-scale neighborhood window set is aligned with the horizontal coordinate of the point to be processed. For each aligned square window, all point cloud data within its boundary are collected to form a neighborhood point set. The arithmetic mean of the elevation values ​​of all points in the neighborhood point set and the standard deviation of these elevation values ​​are calculated. The combination of the arithmetic mean and the standard deviation is recorded as the elevation change feature of the point to be processed at the current window scale. The total number of point clouds in the neighborhood point set is counted, and the total number is divided by the area of ​​the current square window. The result is recorded as the point cloud distribution density feature. The neighborhood point set is projected onto a horizontal plane, and a corresponding two-dimensional grayscale image is generated using the elevation value of the point as the grayscale value. The grayscale co-occurrence matrix of the two-dimensional grayscale image is calculated, and four statistical quantities, energy, contrast, correlation and entropy, are extracted from it. The numerical values ​​of these four statistical quantities are combined and recorded as surface texture features.

5. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 1, characterized in that, In the dynamic filtering module, the process of generating comprehensive feature values ​​is as follows: For a multi-scale visual cognitive feature set of a point within a point cloud data block, traverse each feature value in the set, find the maximum and minimum values ​​of all points on the same feature dimension, and linearly scale the feature value of the current point using the maximum and minimum values ​​to map it to the range of zero to one. Calculate the standard deviation of the elevation values ​​of all points within the point cloud data block, and define the standard deviation as the terrain undulation of the data block. Assign a weight coefficient that is directly proportional to the terrain undulation to the elevation variation feature, assign a weight coefficient that is inversely proportional to the terrain undulation to the surface texture feature, and assign a fixed weight coefficient to the point cloud distribution density feature. Multiply each normalized feature value of the current point by its corresponding weight coefficient and sum them up. The sum is the comprehensive feature value of the point.

6. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 5, characterized in that, In the dynamic filtering module, the process of generating the dynamic threshold interval is as follows: The comprehensive feature values ​​of all points within the point cloud data block are statistically analyzed, and the arithmetic mean and standard deviation of these comprehensive feature values ​​are calculated. The point cloud data block is uniformly divided into multiple rectangular sub-regions on a horizontal plane, and the point cloud density in each sub-region is calculated. The variance of the point cloud density in all sub-regions is calculated, and the variance is defined as the density spatial heterogeneity of the point cloud data block. A preset basic adjustment coefficient is obtained, and the basic adjustment coefficient is multiplied by the density spatial heterogeneity to obtain the dynamic adjustment coefficient. The product of the standard deviation and the dynamic adjustment coefficient is subtracted from the arithmetic mean of the comprehensive feature values, and the result is used as the lower limit of the dynamic threshold interval. The product of the standard deviation and the dynamic adjustment coefficient is added to the arithmetic mean, and the result is used as the upper limit of the dynamic threshold interval.

7. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 1, characterized in that, In the region growth and aggregation module, the process of calculating multi-scale feature similarity is as follows: Obtain the multi-scale visual cognitive feature set for each point in the candidate ground point set. For any two points, calculate the sum of the absolute values ​​of the differences between their corresponding feature values ​​as the feature difference degree. For each point in the set, calculate its feature difference degree with all other points in the set and sum them to obtain the overall difference degree of that point. Sort the candidate ground points in order of increasing overall difference degree. Select the first few points from the sorted sequence to form the initial growth starting point set. For each point in the initial growth starting point set, calculate its feature difference degree with all other points in the candidate ground point set. This feature difference degree is used as the inverse vector degree of the multi-scale feature similarity between the two points.

8. A lidar ground point extraction system based on visual cognitive analog filtering according to claim 7, characterized in that, In the region growth and aggregation module, the process of iteratively growing and aggregating based on the similarity judgment results to form an optimized set of ground points is as follows: Create a global set of processed point labels. Select an unlabeled point from the initial growth starting point set as the current growth seed point. Create an empty current growth subset and add the current growth seed point to it. Obtain the feature difference between the current growth seed point and all unlabeled candidate ground points. Add points with feature difference less than the preset growth similarity threshold to the current growth subset and record these new points as newly added growth seed points in this round. Mark all points in the current growth subset as processed. For each newly added seed point, the process of obtaining feature difference, determining the inclusion condition, and updating the subset and label is repeated until no new points are added. At this point, the current growth subset constitutes an optimized ground point subset. The process of selecting seed points from the initial growth starting point set, creating and growing subsets is repeated until all starting points are marked as processed. Finally, the union of all subsets constitutes the optimized ground point set.

9. The visual cognitive simulation filtering based laser radar ground point extraction system according to claim 1, characterized in that, In the boundary fusion output module, the process of outputting the final ground point extraction result is as follows: Merge the optimized ground point subsets of all point cloud data blocks to form an initial global point cloud. Set a boundary buffer distance and identify all points whose horizontal projected coordinates are less than the boundary buffer distance of any point cloud data block boundary line. Mark them as the boundary point set to be verified. For each point in the boundary point set to be verified, find its neighboring points from different point cloud data blocks and calculate the Euclidean distance between the feature set of the point and the average feature set of the neighboring points. The Euclidean distance is compared with a preset consistency threshold. If the threshold is exceeded, the original feature set of the point is replaced with the average feature value of the neighboring points. For all boundary points whose feature sets have been replaced, their three-dimensional coordinates are recalculated by distance-weighted interpolation based on their new feature values ​​and the coordinates of the neighboring points. After all boundary points are verified and their coordinates are adjusted, the final ground point extraction result is output.

Images