Global localization method for heterogeneous point cloud based on columnar voxel and adaptive down-sampling
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- SOUTHEAST UNIV
- Filing Date
- 2026-04-17
- Publication Date
- 2026-08-07
AI Technical Summary
基于虚拟LiDAR的方法虽能缓解线数不同的差异，但在结构相似场景下，场景表征受限，定位结果不够鲁棒；基于圆柱坐标系径向-高程平面投影的多扇区描述符，虽实现了优异的旋转不变性，但对平移变化敏感；基于自适应网格投影和高斯核加权的方法，生成激光雷达模态不变的表示，但该方法的定位单元为点云子地图，需要累积多帧点云，一次全局定位性能可能受限
[0034]本发明所述的一种基于柱状体素和自适应降采样的异构点云全局定位方法，通过自适应降采样模块自适应确定体素尺寸，同时保留点云几何分布特征；基于柱状体素的特征提取模块，保证特征提取质量的基础上，提高了计算效率，并采用基于拉东变换的描述符提取与位姿估计方法，提高了激光雷达线数异构条件下全局定位的精度、鲁棒性和效率。在异构激光雷达定位的HeLiPR数据集上开展了实验验证，结果表明：所提方法在多个线数异构场景下，相比RING++方法定位成功率有所提升，平均位置精度均保持亚米级，平均姿态精度略高于1°，同时定位效率提升36.3%；为异构多机的自主建图与导航定位提供了可行技术方案。
Smart Images
Figure CN122530296A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of lidar positioning technology, specifically relating to a global positioning method for heterogeneous point clouds based on column voxels and adaptive downsampling. Background Technology
[0002] Localization is a crucial aspect of robotics and autonomous driving, providing prior information for planning and control, and playing a vital role in smart mines, intelligent factories, autonomous parking, and unmanned inspections. However, with the development of LiDAR technology, LiDARs with different line counts, field of view angles, and scanning methods can adapt to various application tasks. But this also means that during localization, maps and queries may come from different types of LiDARs. Traditional LiDAR localization methods are prone to sparse correspondences, ambiguous location recognition, and even localization failure in such scenarios.
[0003] Existing LiDAR localization methods have limitations when dealing with heterogeneous LiDAR line counts. While virtual LiDAR-based methods can mitigate the differences in line counts, they suffer from limited scene representation and insufficient robustness in structurally similar scenarios. Multi-sector descriptors based on radial-elevation plane projection in cylindrical coordinates achieve excellent rotation invariance but are sensitive to translational changes. Methods based on adaptive mesh projection and Gaussian kernel weighting generate mode-invariant representations of the LiDAR, but their localization unit is a point cloud sub-map, requiring the accumulation of multiple point cloud frames, which may limit the performance of a single global localization operation. Therefore, this paper proposes a heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling to address the aforementioned problems. Summary of the Invention
[0004] To address the aforementioned issues, this invention discloses a global localization method for heterogeneous point clouds based on columnar voxels and adaptive downsampling. The method adaptively determines voxel sizes through an adaptive downsampling module while preserving the geometric distribution features of the point cloud. A columnar voxel-based feature extraction module improves computational efficiency while ensuring feature extraction quality. Furthermore, a descriptor extraction and pose estimation method based on Radon transform is employed to enhance the accuracy, robustness, and efficiency of global localization under heterogeneous lidar line count conditions.
[0005] To achieve the above objectives, the technical solution of the present invention is as follows:
[0006] A heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling includes:
[0007] Step 1: The point cloud is segmented into uniformly sized voxels through height and distance filtering, and the feature value of each voxel is calculated. The voxels are then classified based on these feature values. Depending on the current voxel size and classification results, it is determined whether further segmentation of the voxel is needed, until all voxels no longer require segmentation, thus achieving adaptive voxel size determination. Uniform voxel downsampling is then performed on each segmented voxel to achieve adaptive downsampling.
[0008] Step 2: After sampling the point cloud by height, each layer is divided into several cylinders. Information such as the cylinder center and centroid is used to expand the data dimension of the input point cloud. A multilayer perceptron is used to extract features within each cylinder to generate pseudo-images. The pseudo-images from each layer are combined to form a 6-channel BEV feature map.
[0009] Step 3: Using Radon transform, the BEV feature map is converted into a Radon sine map, and then a discrete Fourier transform is used to obtain translation-invariant and rotation-equivariant descriptors. The corresponding descriptors are retrieved from a pre-built map database to achieve location identification. The phase difference between descriptors is calculated using one-dimensional cyclic cross-correlation as a relative heading estimate. Based on rotation compensation, the two-dimensional cross-correlation of the original point cloud is calculated, and its maximum value is used as a relative translation estimate to obtain the global pose.
[0010] The specific steps are as follows:
[0011] Step 1. The point cloud is segmented into uniformly sized voxels through height and distance filtering, and the feature value of each voxel is calculated. The voxels are then classified based on these feature values. Depending on the current voxel size and classification results, it is determined whether further segmentation of the voxel is needed, until all voxels no longer require segmentation, thus achieving adaptive determination of voxel size. Uniform voxel downsampling is then performed on each segmented voxel; specifically as follows:
[0012] The adaptive downsampling method specifically includes three steps: First, the point cloud is downsampled according to a set side length. The data is divided into several uniformly sized voxels. The covariance matrix of each voxel is calculated. Then, eigenvalue decomposition is performed to obtain three eigenvalues, and these eigenvalues are defined. Based on the eigenvalues, the three geometric model factors are calculated according to the following formula:
[0013] (1);
[0014] Where represents the linear model factor, and represents the planar model factor. This represents the solid model factor. By comparing the maximum values of the three geometric model factors, the voxels are divided into three geometric models: linear, planar, and solid models, i.e.:
[0015] (2);
[0016] when When =1, it means The value of is the largest, therefore the voxel is classified as a linear voxel model. Similarly, when When the value is 2 and 3, it means that the voxel is classified as a planar model and a solid model, respectively.
[0017] The characteristic of the linear voxel model is that the point cloud is mainly distributed along a specific direction, and this direction can be determined by the covariance matrix of the voxels. The process is as follows: Project all points along the main direction. By measuring the distance between each projected point and its two adjacent points, the point cloud distribution density along that direction can be evaluated. If there is uneven point cloud distribution, the spaces separating the dense regions will be considered as sparse regions within the voxel. The centroid c is determined as the dividing point, and dividing lines are arranged along the x, y, and z axes and pass through this point, dividing the voxel into 8 blocks, thus achieving adaptive determination of the voxel size.
[0018] To avoid incorrect segmentation, the distance between the points and the dividing line is measured before segmentation. If the proportion of points close to the dividing line exceeds a threshold... If uniform downsampling is used, the result is acceptable, and therefore the segmentation result will be discarded. Adaptive voxel segmentation followed by uniform voxel downsampling reduces the number of points while preserving their geometric distribution characteristics.
[0019] The segmentation operation for planar models is similar to that for linear models, the main difference being that voxels in planar models are projected onto a plane defined by two principal directions. Sparse regions are then extracted from the projection plane to complete voxel segmentation. For 3D models, considering their complex spatial distribution, the voxel is divided into eight equal-sized sub-voxels and reclassified until the voxel radius is less than a threshold. Therefore, this invention only considers linear voxel models and planar voxel models. These three steps are repeated until all voxels no longer require segmentation.
[0020] Step 2. The point cloud sampled according to height in this invention is divided into several cylinders in each layer. Information such as the cylinder center and centroid is combined to expand the data dimension of the input point cloud. A multilayer perceptron is used to extract features within each cylinder to generate pseudo-images. The pseudo-images from each layer of the point cloud are combined to form a 6-channel BEV feature map; specifically as follows:
[0021] The input point cloud of the feature extraction module is denoted as To further improve the efficiency of BEV extraction, this invention improves the Pillar Feature Net (PFN) and Voxel Feature Encoding (VFE) modules in PointPillars, achieving efficient and reliable feature extraction. The PFN module divides the point cloud into a voxel grid composed of several groups of voxels, and augments the original points using information such as the centroid of the pillars. The VFE module maps the augmented point cloud features to a high-dimensional feature space, then selects appropriate features to convert them into BEVs and outputs them.
[0022] Before feature extraction, the downsampled point cloud is first sorted according to a distance threshold. ]、[ ] and height threshold[ Distant points and ground points are removed separately, and the point cloud is then divided into 6 layers along the height dimension, with each layer corresponding to a feature extraction channel. For each layer, the PFN module divides the point cloud along the xz and yz planes according to the side lengths. The data is divided into multiple cylinders. Each cylinder has a square cross-section and extends vertically to the boundary. Since the original coordinate distribution of the point cloud is related to the observation location, which affects the translation invariance of subsequent descriptor extraction, and the reflectivity intensity information characteristics may differ between different LiDAR models, the original coordinates and reflectivity intensity information of the point cloud within each non-empty cylinder are not included in feature extraction. Instead, the input features of each point are expanded to [xc, yc, zc, xp, yp, zp] to achieve data augmentation. Here, [xc, yc, zc] represents the coordinates of the cylinder's centroid, and [xp, yp, zp] represents the offset of the point from the cylinder's center. Due to the sparsity of the point cloud, there are a large number of empty cylinders, affecting feature extraction efficiency. Therefore, the sparse cylinders are reorganized into a tensor of size [D, P, N], where D represents the feature dimension of each point, P represents the number of non-empty cylinders, and N represents the number of points in each cylinder. For a non-empty pillar, if the number of point clouds within the pillar exceeds N, then the point clouds are randomly sampled; if the number of point clouds is less than N, then zeros are used to pad the points.
[0023] The feature encoding module inputs the tensor obtained from feature extraction into the feature network, and uses a multilayer perceptron network to encode the enhanced data, expanding the original D-dimensional features to C-dimensional.
[0024] First, a linear layer is applied to each point, followed by batch normalization and rectified linear unit (ReLU) to obtain a tensor of [C, P, N] dimensions:
[0025] (3);
[0026] In the formula, Let represent a three-dimensional tensor of size [D, P, N]. Represents the ReLU activation function. This indicates a batch normalization operation. This represents a three-dimensional tensor of size [C, P, N] obtained after passing through a multilayer perceptron.
[0027] Subsequently, the maximum value of the tensor within each cylinder is taken as the eigenvalue of that cylinder, and the eigenvalues of all cylinders are placed back into the coordinates of the cylinder to obtain the BEV representation of the point cloud of that layer, denoted as . The BEV feature map obtained by combining the BEVs of each layer is denoted as . This serves as the global feature representation of the current point cloud frame. The method of this invention improves the efficiency of feature extraction through dense tensor representation, while the fusion of highly hierarchical extraction preserves the point cloud distribution information along the z-axis.
[0028] Step 3. Using Radon transform, the BEV feature map is converted into a Radon sine map, and then a translation-invariant and rotation-equivariant descriptor is obtained through discrete Fourier transform. The corresponding descriptor is retrieved from a pre-built map database to achieve location identification. The phase difference between descriptors is calculated using one-dimensional cyclic cross-correlation as a relative heading estimate. Based on rotation compensation, the two-dimensional cross-correlation of the original point cloud is calculated, and its maximum value is used as a relative translation estimate to obtain the global pose. Details are as follows:
[0029] After completing the BEV extraction, this invention applies the Radon transform to the BEV, aggregates it to obtain a globally consistent descriptor with rotation-translation invariance, and estimates the 3-DOF global pose based on the descriptor and BEV.
[0030] First, for each layer of BEV Perform Radon transforms on each part to obtain the Radon sine curve. And the Radon sine diagram along The shaft exhibits rotational isotropy, and along... non-uniform shaft displacement To completely eliminate non-uniform displacement, a row-by-row one-dimensional DFT is used along the Radon sine curve. The axis is transformed to the frequency domain, and then the amplitude of the Fourier spectrum is calculated. The resulting representation has translation invariance and is called a translation-invariant-rotation equivariant graph. Variant Rotation-Equivariant Gram (TING), denoted as TING. Due to the decoupling of the transformation dimensions, DFT eliminates the changes caused by translation and does not violate rotational equivariance; therefore, TING possesses translation invariance and rotational equivariance. The 6-channel BEV is transformed using the above method and then concatenated to obtain the 6-channel TING.
[0031] For query point cloud Q and map point cloud M i Extract the descriptor TING separately, and denote TING for the query point cloud as and TING for the map point cloud as . First, TING representation of each map point cloud Mi in the map database ,along By performing a one-dimensional cyclic cross-correlation along the axis and accumulating the results along the channel dimension, a single-channel one-dimensional cross-correlation map can be obtained, ultimately generating a batch of cross-correlation maps. Due to TING's translation invariance, The axes have no effect in the two-dimensional cross-correlation calculation. Therefore, the two-dimensional cross-correlation between two TINGs can be simplified to along... One-dimensional cross-correlation of axes. Subsequently, for each TING pair ( , Related images Max pooling is performed, and the resulting descriptor is denoted as... Due to the invariance of the maximum value with respect to attitude change, and Since they are the same, RING has rotation-translation invariance.
[0032] The pose estimation stage comprises four steps: position identification, rotation estimation, translation estimation, and refinement. This invention first uses a global descriptor, the RING, which is rotation-translation invariant, for position identification. However, due to occlusion effects, the invariance of the RING cannot be precisely guaranteed when the translation amount changes significantly relative to the scan range. Therefore, by comparing the RING differences between the query scan and the map scan, the closest map scan can be retrieved from the map database. This enables location identification. After location identification is completed, the query scan is accumulated. With retrieved map scan Correlation diagram between them This allows us to obtain an estimate of the relative rotation. Based on this estimate of the relative rotation, to further estimate the relative translation, this invention first... Rotation - To compensate for relative rotation, the correlation diagram for each channel is then calculated. Then, the correlation maps are accumulated along the channel dimension, and the maximum value of the correlation map is used as the estimate of the global translation. In the pose refinement part, the estimated value of the 3-DOF relative transformation between the query point cloud and the map point cloud is used as the initial parameters of the ICP algorithm, and finally the optimized 6-DOF pose is obtained.
[0033] The beneficial effects of this invention are as follows:
[0034] This invention presents a heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling. The method adaptively determines voxel sizes through an adaptive downsampling module while preserving the geometric distribution features of the point cloud. A columnar voxel-based feature extraction module improves computational efficiency while ensuring feature extraction quality. Furthermore, a descriptor extraction and pose estimation method based on Radon transform enhances the accuracy, robustness, and efficiency of global localization under heterogeneous LiDAR line count conditions. Experimental verification was conducted on the HeLiPR dataset for heterogeneous LiDAR localization. Results show that the proposed method improves the localization success rate compared to the RING++ method in multiple heterogeneous line count scenarios, maintaining sub-meter average position accuracy and slightly exceeding 1° average pose accuracy, while improving localization efficiency by 36.3%. This provides a feasible technical solution for autonomous mapping and navigation localization in heterogeneous multi-LiDAR environments. Attached Figure Description
[0035] Figure 1 This is a flowchart illustrating the implementation of the heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling as described in this invention.
[0036] Figure 2 This is a comparison chart of the pose accuracy of the method described in this invention compared to other open-source methods under heterogeneous conditions;
[0037] Figure 3 This is a time consumption statistics chart of the method described in this invention. Detailed Implementation
[0038] The present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments. It should be understood that the following specific embodiments are for illustrative purposes only and are not intended to limit the scope of the invention.
[0039] like Figure 1 As shown, the heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling described in this invention comprises the following specific steps:
[0040] Step 1: The point cloud is segmented into uniformly sized voxels through height and distance filtering, and the feature value of each voxel is calculated. The voxels are then classified based on these feature values. Depending on the current voxel size and classification results, it is determined whether further segmentation of the voxel is needed, until all voxels no longer require segmentation, thus achieving adaptive determination of voxel size. Uniform voxel downsampling is then performed on each segmented voxel.
[0041] Specifically as follows:
[0042] The adaptive downsampling method specifically includes three steps: First, the point cloud is downsampled according to a set side length. The data is divided into several uniformly sized voxels. The covariance matrix of each voxel is calculated. Then, eigenvalue decomposition is performed to obtain eigenvalues, and eigenvalues are defined. Based on the eigenvalues, the three geometric model factors are calculated according to the following formula:
[0043] (1);
[0044] in, Represents the linear model factor. Represents the plane model factor. This represents the solid model factor. By comparing the maximum values of the three geometric model factors, the voxels are divided into three geometric models: linear, planar, and solid models, i.e.:
[0045] (2);
[0046] when When =1, it means The value of is the largest, therefore the voxel is classified as a linear voxel model. Similarly, when When the value is 2 and 3, it means that the voxel is classified as a planar model and a solid model, respectively.
[0047] The characteristic of the linear voxel model is that the point cloud is mainly distributed along a specific direction, and this direction can be determined by the covariance matrix of the voxels. The process is as follows: Project all points along the main direction. By measuring the distance between each projected point and its two adjacent points, the point cloud distribution density along that direction can be evaluated. If there is uneven point cloud distribution, the spaces separating the dense regions will be considered as sparse regions within the voxel. The centroid c is determined as the dividing point, and dividing lines are arranged along the x, y, and z axes and pass through this point, dividing the voxel into 8 blocks, thus achieving adaptive determination of the voxel size.
[0048] To avoid incorrect segmentation, the distance between the points and the dividing line is measured before segmentation. If the proportion of points close to the dividing line exceeds a threshold... If uniform downsampling is used, the result is acceptable, and therefore the segmentation result will be discarded. Adaptive voxel segmentation followed by uniform voxel downsampling reduces the number of points while preserving their geometric distribution characteristics.
[0049] The segmentation operation for planar models is similar to that for linear models, the main difference being that voxels in planar models are projected onto a plane defined by two principal directions. Sparse regions are then extracted from the projection plane to complete voxel segmentation. For 3D models, considering their complex spatial distribution, the voxel is divided into eight equal-sized sub-voxels and reclassified until the voxel radius is less than a threshold. Therefore, this invention only considers linear voxel models and planar voxel models. These three steps are repeated until all voxels no longer require segmentation.
[0050] Step 2: After sampling the point cloud by height, each layer is divided into several cylinders. Information such as the cylinder center and centroid is used to expand the data dimension of the input point cloud. A multilayer perceptron is used to extract features within each cylinder to generate pseudo-images. The pseudo-images from each layer are combined to form a 6-channel BEV feature map.
[0051] Specifically as follows:
[0052] The input point cloud of the feature extraction module is denoted as To further improve the efficiency of BEV extraction, this invention improves the Pillar Feature Net (PFN) and Voxel Feature Encoding (VFE) modules in PointPillars, achieving efficient and reliable feature extraction. The PFN module divides the point cloud into a voxel grid composed of several groups of voxels, and augments the original points using information such as the centroid of the pillars. The VFE module maps the augmented point cloud features to a high-dimensional feature space, then selects appropriate features to convert them into BEVs and outputs them.
[0053] Before feature extraction, the downsampled point cloud is first sorted according to a distance threshold. ]、[ ] and height threshold[ Distant points and ground points are removed separately, and the point cloud is then divided into 6 layers along the height dimension, with each layer corresponding to a feature extraction channel. For each layer, the PFN module divides the point cloud along the xz and yz planes according to the side lengths. The data is divided into multiple cylinders. Each cylinder has a square cross-section and extends vertically to the boundary. Since the original coordinate distribution of the point cloud is related to the observation location, which affects the translation invariance of subsequent descriptor extraction, and the reflectivity intensity information characteristics may differ between different LiDAR models, the original coordinates and reflectivity intensity information of the point cloud within each non-empty cylinder are not included in feature extraction. Instead, the input features of each point are expanded to [xc, yc, zc, xp, yp, zp] to achieve data augmentation. Here, [xc, yc, zc] represents the coordinates of the cylinder's centroid, and [xp, yp, zp] represents the offset of the point from the cylinder's center. Due to the sparsity of the point cloud, there are a large number of empty cylinders, affecting feature extraction efficiency. Therefore, the sparse cylinders are reorganized into a tensor of size [D, P, N], where D represents the feature dimension of each point, P represents the number of non-empty cylinders, and N represents the number of points in each cylinder. For a non-empty pillar, if the number of point clouds within the pillar exceeds N, then the point clouds are randomly sampled; if the number of point clouds is less than N, then zeros are used to pad the points.
[0054] The feature encoding module inputs the tensor obtained from feature extraction into the feature network, and uses a multilayer perceptron network to encode the enhanced data, expanding the original D-dimensional features to C-dimensional.
[0055] First, a linear layer is applied to each point, followed by batch normalization and rectified linear unit (ReLU) to obtain a tensor of [C, P, N] dimensions:
[0056] (3);
[0057] In the formula, Let represent a three-dimensional tensor of size [D, P, N]. Represents the ReLU activation function. This indicates a batch normalization operation. This represents a three-dimensional tensor of size [C, P, N] obtained after passing through a multilayer perceptron.
[0058] Subsequently, the maximum value of the tensor within each cylinder is taken as the eigenvalue of that cylinder, and the eigenvalues of all cylinders are placed back into the coordinates of the cylinder to obtain the BEV representation of the point cloud of that layer, denoted as . The BEV feature map obtained by combining the BEVs of each layer is denoted as . This serves as the global feature representation of the current point cloud frame. The method of this invention improves the efficiency of feature extraction through dense tensor representation, while the fusion of highly hierarchical extraction preserves the point cloud distribution information along the z-axis.
[0059] Step 3: Using Radon transform, the BEV feature map is converted into a Radon sine map, and then a discrete Fourier transform is used to obtain translation-invariant and rotation-equivariant descriptors. The corresponding descriptors are retrieved from a pre-built map database to achieve location identification. The phase difference between descriptors is calculated using one-dimensional cyclic cross-correlation as a relative heading estimate. Based on rotation compensation, the two-dimensional cross-correlation of the original point cloud is calculated, and its maximum value is used as a relative translation estimate to obtain the global pose.
[0060] Specifically as follows:
[0061] After completing the BEV extraction, this invention applies the Radon transform to the BEV, aggregates it to obtain a globally consistent descriptor with rotation-translation invariance, and estimates the 3-DOF global pose based on the descriptor and BEV.
[0062] First, for each layer of BEV Perform Radon transforms on each part to obtain the Radon sine curve. And the Radon sine diagram along The shaft exhibits rotational isotropy, and along... non-uniform shaft displacement To completely eliminate non-uniform displacement, a row-by-row one-dimensional DFT is used along the Radon sine curve. The axis is transformed to the frequency domain, and then the amplitude of the Fourier spectrum is calculated. The resulting representation has translation invariance and is called the Translation-Invariant Rotation-Equivariant Gram (TING), denoted as . Due to the decoupling of the transformation dimensions, the DFT eliminates the changes caused by translation and does not violate rotational isotropy; therefore, TING possesses translation invariance and rotational isotropy. The 6-channel BEV is transformed using the above methods and then concatenated to obtain the 6-channel TING.
[0063] For query point cloud Q and map point cloud M i Extract the descriptor TING separately, and denote the query point cloud TING as The TING of the map point cloud is First, TING representation of each map point cloud Mi in the map database ,along By performing a one-dimensional cyclic cross-correlation along the axis and accumulating the results along the channel dimension, a single-channel one-dimensional cross-correlation map can be obtained, ultimately generating a batch of cross-correlation maps. Due to TING's translation invariance, The axes have no effect in the two-dimensional cross-correlation calculation. Therefore, the two-dimensional cross-correlation between two TINGs can be simplified to along... One-dimensional cross-correlation of axes. Subsequently, for each TING pair ( , Related images Max pooling is performed, and the resulting descriptor is denoted as... Due to the invariance of the maximum value with respect to attitude change, and Since they are the same, RING has rotation-translation invariance.
[0064] The pose estimation stage comprises four steps: position identification, rotation estimation, translation estimation, and refinement. This invention first uses a global descriptor, the RING, which is rotation-translation invariant, for position identification. However, due to occlusion effects, the invariance of the RING cannot be precisely guaranteed when the translation amount changes significantly relative to the scan range. Therefore, by comparing the RING differences between the query scan and the map scan, the closest map scan can be retrieved from the map database. This enables location identification. After location identification is completed, the query scan is accumulated. With retrieved map scan Correlation diagram between them This allows us to obtain an estimate of the relative rotation. Based on this estimate of the relative rotation, to further estimate the relative translation, this invention first... Rotation - To compensate for relative rotation, the correlation diagram for each channel is then calculated. Then, the correlation maps are accumulated along the channel dimension, and the maximum value of the correlation map is used as the estimate of the global translation. In the pose refinement part, the estimated value of the 3-DOF relative transformation between the query point cloud and the map point cloud is used as the initial parameters of the ICP algorithm, and finally the optimized 6-DOF pose is obtained.
[0065] It should be noted that the above content merely illustrates the technical concept of the present invention and should not be construed as limiting the scope of protection of the present invention. For those skilled in the art, various improvements and modifications can be made without departing from the principle of the present invention, and all such improvements and modifications fall within the scope of protection of the claims of the present invention.
Claims
1. A heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling, characterized in that, include: Step 1: The point cloud is segmented into uniformly sized voxels by applying height and distance filtering, and the feature value of each voxel is calculated. The voxels are then classified based on the feature value. Based on the current voxel size and classification results, it is determined whether to continue segmenting the voxel until all voxels no longer need to be segmented, thus achieving adaptive determination of voxel size. Uniform voxel downsampling is then performed on each segmented voxel. Step 2: After sampling the point cloud by height, each layer of the point cloud is divided into several pillars. The data dimension of the input point cloud is expanded by combining the information of the pillar center and centroid. Features are extracted in each pillar using a multilayer perceptron to generate pseudo-images. The pseudo-images of each layer of the point cloud are combined into a 6-channel BEV feature map. Step 3: Using Radon transform, the BEV feature map is converted into a Radon sine map, and then a translation-invariant and rotation-equivariant descriptor is obtained through discrete Fourier transform; the corresponding descriptor is retrieved from the pre-built map database to achieve position recognition; the phase difference between the descriptors is calculated by one-dimensional cyclic cross-correlation as a relative heading estimate; based on the compensated rotation, the two-dimensional cross-correlation of the original point cloud is calculated, and its maximum value is used as a relative translation estimate to obtain the global pose.
2. The heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling as described in claim 1, characterized in that, Step 1 involves applying height and distance filters to segment the point cloud into uniformly sized voxels, calculating the feature value for each voxel, and then classifying the voxels based on these feature values. Based on the current voxel size and classification results, a decision is made as to whether to continue segmenting the voxel until all voxels no longer require segmentation, thus achieving adaptive determination of voxel size. Uniform voxel downsampling is then performed on each segmented voxel. The specific steps are as follows: The adaptive downsampling method specifically includes three steps: First, the point cloud is downsampled according to a set side length. The data is divided into several uniformly sized voxels; the covariance matrix of each voxel is calculated. Then, eigenvalue decomposition is performed to obtain eigenvalues, and eigenvalues are defined. Based on the eigenvalues, the three geometric model factors are calculated according to the following formula: （1）； in, Represents the linear model factor. Represents the plane model factor. This represents the solid model factor; by comparing the maximum values of the three geometric model factors, the voxels are divided into three geometric models: linear, planar, and solid models, i.e.: (2)； when When =1, it means The value of is the largest, therefore the voxel is classified as a linear voxel model; similarly, when When the value is 2 and 3, it means that the voxel is classified as a planar model and a solid model, respectively. To avoid incorrect segmentation, the distance between points and the dividing line is measured before segmentation; if the proportion of points close to the dividing line exceeds a threshold... If so, discard the segmentation result; The segmentation operation for planar models is similar to that for linear models. Voxels are projected onto a plane defined by two principal directions. Sparse regions are then extracted from the projection plane to complete the voxel segmentation. For 3D models, voxels are divided into 8 equal-sized sub-voxels and then reclassified until the voxel radius is less than a threshold. The above three steps are repeated until all voxels no longer need to be segmented.
3. The heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling as described in claim 1, characterized in that, Step 2 involves sampling the point cloud by height, with each layer divided into several cylinders. The data dimensions of the input point cloud are expanded by combining the cylinder center and centroid information. Features are extracted within each cylinder using a multilayer perceptron to generate pseudo-images. These pseudo-images from each layer are then combined to form a 6-channel BEV feature map. The details are as follows: The input point cloud of the feature extraction module is denoted as ; The columnar voxel feature extraction and voxel feature encoding modules in PointPillars have been improved. The PFN module divides the point cloud into a voxel grid composed of several groups of voxels and uses the centroid information of the columns to augment the original points. The VFE module maps the augmented point cloud features to a high-dimensional feature space, and then selects appropriate features to convert them into BEV and output them. Before feature extraction, the downsampled point cloud is first sorted according to a distance threshold. ]、[ ] and height threshold[ [Remove distant points and ground points separately, then divide the point cloud into 6 layers along the height dimension, with each layer corresponding to a feature extraction channel; for each layer, the PFN module divides the point cloud along the xz and yz planes according to the side length] The data is divided into multiple cylindrical sections with square cross-sections and vertical heights extending to the boundaries. The input features of each point are then expanded to [xc, yc, zc, xp, yp, zp] to augment the data. Here, [xc, yc, zc] represents the coordinates of the cylindrical centroid, and [xp, yp, zp] represents the offset of the point from the cylindrical centroid. The sparse cylindrical sections are then reorganized into tensors of size [D, P, N], where D represents the feature dimension of each point, P represents the number of non-empty cylindrical sections, and N represents the number of points in each cylindrical section. For non-empty cylindrical sections, if the number of point clouds within a section exceeds N, the point clouds are randomly sampled; if the number of point clouds is less than N, zeros are used for padding. The feature encoding module inputs the tensor obtained from feature extraction into the feature network, and uses a multilayer perceptron network to encode the enhanced data, expanding the original D-dimensional features to C-dimensional; First, apply a linear layer to each point, and then sequentially pass through batch normalization and corrected linear units to obtain a tensor of [C, P, N] dimensions: (3)； In the formula, Let represent a three-dimensional tensor of size [D, P, N]. Represents the ReLU activation function. This indicates a batch normalization operation. This represents a three-dimensional tensor of size [C, P, N] obtained after passing through a multilayer perceptron. Subsequently, the maximum value of the tensor within each cylinder is taken as the eigenvalue of that cylinder, and the eigenvalues of all cylinders are placed back into the coordinates of the cylinder to obtain the BEV representation of the point cloud of that layer, denoted as . The BEV feature map obtained by combining the BEVs of each layer is denoted as . This serves as the global feature representation of the current point cloud frame.
4. The heterogeneous point cloud global localization method based on columnar voxels and adaptive downsampling as described in claim 1, characterized in that, Using Radon transform, the BEV feature map is converted into a Radon sine map, and then a translation-invariant and rotation-equivalent descriptor is obtained through discrete Fourier transform. The corresponding descriptor is retrieved from a pre-built map database to achieve location identification. The phase difference between the descriptors is calculated using one-dimensional cyclic cross-correlation as a relative heading estimate. Based on rotation compensation, the two-dimensional cross-correlation of the original point cloud is calculated, and its maximum value is used as a relative translation estimate to obtain the global pose. The details are as follows: After BEV extraction, Radon transform is applied to BEV to aggregate globally consistent descriptors with rotation-translation invariance. Based on the descriptors and BEV, the 3-DOF global pose is estimated. First, for each layer of BEV Perform Radon transforms on each part to obtain the Radon sine curve. And the Radon sine diagram along The shaft exhibits rotational isotropy, and along... non-uniform shaft displacement To completely eliminate non-uniform displacement, a row-by-row one-dimensional DFT is used along the Radon sine curve. The axis is transformed to the frequency domain, and then the amplitude of the Fourier spectrum is calculated. The resulting representation has translation invariance and is called a translation-invariant-rotational isomorphic graph, denoted as . The 6-channel BEV is transformed as described above and then spliced together to obtain the 6-channel TING. For query point cloud Q and map point cloud M i Extract the descriptor TING separately, and denote the query point cloud TING as The TING of the map point cloud is First, TING representation of each map point cloud Mi in the map database ,along One-dimensional cyclic cross-correlation is performed along the axis, and the results are accumulated along the channel dimension to obtain a one-dimensional cross-correlation map for a single channel, ultimately generating a batch of cross-correlation maps; due to the translation invariance of TING, The axes have no effect in two-dimensional cross-correlation calculations; Therefore, the two-dimensional cross-correlation between the two TINGs simplifies to along One-dimensional cross-correlation of axes; Subsequently, for each TING pair ( , Related images Max pooling is performed, and the resulting descriptor is denoted as... Due to the invariance of the maximum value with respect to attitude change, and The same applies, therefore RING has rotation-translation invariance; The pose estimation stage comprises four steps: position identification, rotation estimation, translation estimation, and refinement. First, position identification is performed using a global descriptor RING, which is rotation-translation invariant. Then, by comparing the RING differences between the query scan and the map scan, the closest map scan is retrieved from the map database. This enables location identification; after location identification is completed, the query scan is accumulated. With retrieved map scan Correlation diagram between them This yields an estimate of the relative rotation; based on the estimated relative rotation, to further estimate the relative translation, we first... Rotation - To compensate for relative rotation, the correlation diagram for each channel is then calculated. Then, the correlation graphs are accumulated along the channel dimension, and the maximum value of the correlation graphs is used as the estimate of the global translation. In the pose refinement part, the estimated value of the 3-DOF relative transformation between the query point cloud and the map point cloud is used as the initial parameter of the ICP algorithm, and finally the optimized 6-DOF pose is obtained.