点云降采样
第一章 点云数据采集
第二章 点云滤波
第二章 点云降采样
1. 为什么要降采样?
我们获得的数据量大,特别是几十万个以上的点云,里面有很多冗余数据,会导致处理起来比较耗时。
降采样是一种有效的减少数据、缩减计算量的方法。
2.降采样算法
2.1 随机降采样
根据设置的比例系数随机删除点云,比较接近均匀采样,但不稳定。
Open3d
import numpy as np
import open3d as o3dpcd = o3d.io.read_point_cloud("second_radius_cloud.pcd")
print(pcd) # 输出点云点的个数
o3d.visualization.draw_geometries([pcd], window_name="原始点云",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
downpcd = pcd.random_down_sample(sampling_ratio=0.5)
print(downpcd) #降采样后的点云数
o3d.visualization.draw_geometries([downpcd], window_name="随机降采样",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
PCL
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/filters/random_sample.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("1697165371469.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;pcl::RandomSample<pcl::PointXYZ> random_sampling;random_sampling.setInputCloud(cloud);random_sampling.setSample(10000); // 设置希望得到的点数random_sampling.filter(*cloud_downsampled);std::cout << "downsampled cloud size: " << cloud_downsampled->width * cloud_downsampled->height << std::endl;pcl::visualization::PCLVisualizer::Ptr viewer(new pcl::visualization::PCLVisualizer("3D Viewer"));viewer->setBackgroundColor(0, 0, 0); // 设置背景色viewer->addPointCloud<pcl::PointXYZ>(cloud_downsampled, "sample cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "sample cloud");viewer->initCameraParameters();while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}
2.2 均匀降采样
就是每隔多远采集一个点,
Open3d
import numpy as np
import open3d as o3dpcd = o3d.io.read_point_cloud("second_radius_cloud.pcd")
print(pcd) # 输出点云点的个数
o3d.visualization.draw_geometries([pcd], window_name="原始点云",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
downpcd = pcd.uniform_down_sample(6)
print(downpcd) #降采样后的点云数
o3d.visualization.draw_geometries([downpcd], window_name="均匀降采样",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
PCL
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/filters/uniform_sampling.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("1697165371469.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;pcl::UniformSampling<pcl::PointXYZ> filter; // 创建均匀采样对象filter.setInputCloud(cloud); // 设置待采样点云filter.setRadiusSearch(10.0f); // 设置采样半径filter.filter(*cloud_downsampled); // 执行均匀采样,结果保存在cloud_filtered中std::cout << "downsampled cloud size: " << cloud_downsampled->width * cloud_downsampled->height << std::endl;pcl::visualization::PCLVisualizer::Ptr viewer(new pcl::visualization::PCLVisualizer("3D Viewer"));viewer->setBackgroundColor(0, 0, 0); // 设置背景色viewer->addPointCloud<pcl::PointXYZ>(cloud_downsampled, "sample cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "sample cloud");viewer->initCameraParameters();viewer->saveScreenshot("screenshot.png");while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}
2.3 体素降采样
将空间切割为均匀大小的体素网格,以非空体素的质心代替该体素内的所有点。
原点云位置使用体素降采样后会发生变化。
open3d
import numpy as np
import open3d as o3d
pcd = o3d.io.read_point_cloud("second_radius_cloud.pcd")
print(pcd) # 输出点云点的个数
o3d.visualization.draw_geometries([pcd], window_name="原始点云",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)downpcd = pcd.voxel_down_sample(voxel_size=5)
print(downpcd)
o3d.visualization.draw_geometries([downpcd], window_name="体素降采样",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
pcl
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/filters/voxel_grid.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("1697165371469.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;pcl::VoxelGrid<pcl::PointXYZ> sor;sor.setInputCloud(cloud);sor.setLeafSize(10.0f, 10.0f, 10.0f);sor.filter(*cloud_downsampled);std::cout << "downsampled cloud size: " << cloud_downsampled->width * cloud_downsampled->height << std::endl;pcl::visualization::PCLVisualizer::Ptr viewer(new pcl::visualization::PCLVisualizer("3D Viewer"));viewer->setBackgroundColor(0, 0, 0); // 设置背景色viewer->addPointCloud<pcl::PointXYZ>(cloud_sampled, "sample cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "sample cloud");viewer->initCameraParameters();viewer->saveScreenshot("screenshot.png");while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}
2.4 最远点降采样
首先随机选择一个点,其次,在剩下点中寻找最远的点,再去再剩下点中找到同时离这两个点最远的点,直到满足采样点个数。
Open3d
import numpy as np
import open3d as o3dpcd = o3d.io.read_point_cloud("second_radius_cloud.pcd")
print(pcd) # 输出点云点的个数
o3d.visualization.draw_geometries([pcd], window_name="原始点云",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
downpcd=pcd.farthest_point_down_sample(10000)
print(downpcd) #降采样后的点云数
o3d.visualization.draw_geometries([downpcd], window_name="最远点降采样",width=1024, height=768,left=50, top=50,mesh_show_back_face=True)
PCL
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/common/distances.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("1697165371469.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;size_t N = cloud->size();assert(N >= 10000);srand(time(0));size_t seed_index = rand() % N;pcl::PointXYZ p = cloud->points[seed_index];;cloud_downsampled->push_back(p);cloud->erase(cloud->begin() + seed_index);for (size_t i = 1; i < 10000; i++){float max_distance = 0;size_t max_index = 0;for (size_t j = 0; j < cloud->size(); j++){float distance = pcl::euclideanDistance(p, cloud->points[j]);if (distance > max_distance){max_distance = distance;max_index = max_index;}}p = cloud->points[max_index];cloud_downsampled->push_back(p);cloud->erase(cloud->begin() + max_index);}std::cout << "downsampled cloud size: " << cloud_downsampled->width * cloud_downsampled->height << std::endl;pcl::visualization::PCLVisualizer::Ptr viewer(new pcl::visualization::PCLVisualizer("3D Viewer"));viewer->setBackgroundColor(0, 0, 0); // 设置背景色viewer->addPointCloud<pcl::PointXYZ>(cloud_downsampled, "sample cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "sample cloud");viewer->initCameraParameters();while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}
2.5 移动最小二乘法降采样
在MLS法中,需要在一组不同位置的节点附近建立拟合曲线,每个节点都有自己的一组系数用于定义该位置附近拟合曲线的形态。因此,在计算某个节点附近的拟合曲线时,只需要计算该点的该组系数值即可。
此外,每个节点的系数取值只考虑其临近采样点,且距离节点越近的采样点贡献越大,对于未置较远的点则不予考虑。
许多文章都将移动最小二乘法作为降采样方法,我觉得这只是一种平滑,所以这里给了重建代码,不进一步实验了。
PCL
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/filters/voxel_grid.h>
#include <pcl/surface/mls.h>
#include <pcl/search/kdtree.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("second_radius_cloud.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;pcl::search::KdTree<pcl::PointXYZ>::Ptr tree (new pcl::search::KdTree<pcl::PointXYZ>);// 输出的PointCloud中有PointNormal类型,用来存储MLS算出的法线pcl::PointCloud<pcl::PointNormal> mls_points;// 定义MovingLeastSquares对象并设置参数pcl::MovingLeastSquares<pcl::PointXYZ, pcl::PointNormal> mls;mls.setComputeNormals(true);mls.setInputCloud(cloud);mls.setSearchMethod(tree);mls.setSearchRadius(30);// 曲面重建mls.process(mls_points);//std::cout << "downsampled cloud size: " << mls_points->width * mls_points->height << std::endl;// 使用PCLVisualizer进行可视化boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer(new pcl::visualization::PCLVisualizer("MLS Cloud Viewer"));viewer->addPointCloud<pcl::PointNormal>(mls_points.makeShared(), "MLS Cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "MLS Cloud");viewer->addPointCloudNormals<pcl::PointNormal>(mls_points.makeShared(), 1, 0.05, "normals"); // 可选:显示法线viewer->saveScreenshot("screenshot.png");while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}
2.6 法线空间采样
通过在法向量空间内均匀随机抽样,使所选点之间的法线分布尽可能大,结果表现为地物特征变化大的地方剩余点较多,变化小的地方剩余点稀少,可有效保持地物特征。
Open3d
import open3d as o3d
import numpy as npdef normal_space_sampling(pcd, num_bins=5, num_samples=10000):# 1. 估算法线pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius=10, max_nn=30))normals = np.asarray(pcd.normals)# 2. 使用法线的x、y和z分量将法线映射到一个3D直方图或“bin”空间bins = np.linspace(-1, 1, num_bins)normal_bins = np.digitize(normals, bins)unique_bins = np.unique(normal_bins, axis=0)sampled_indices = []for b in unique_bins:indices = np.all(normal_bins == b, axis=1)bin_points = np.where(indices)[0]if bin_points.size > 0:sampled_indices.append(np.random.choice(bin_points))# 如果采样点数不足,从原点云中随机选择其他点while len(sampled_indices) < num_samples:sampled_indices.append(np.random.randint(0, len(pcd.points)))# 3. 从每个bin中选择一个点进行采样sampled_points = np.asarray(pcd.points)[sampled_indices]sampled_pcd = o3d.geometry.PointCloud()sampled_pcd.points = o3d.utility.Vector3dVector(sampled_points)return sampled_pcd# 读取点云
pcd = o3d.io.read_point_cloud("second_radius_cloud.pcd")
sampled_pcd = normal_space_sampling(pcd)
o3d.visualization.draw_geometries([sampled_pcd])
PCL
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <pcl/filters/normal_space.h>
#include <pcl/features/normal_3d.h>int main(int argc, char** argv) {pcl::PointCloud<pcl::PointXYZ>::Ptr cloud (new pcl::PointCloud<pcl::PointXYZ>);pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_downsampled(new pcl::PointCloud<pcl::PointXYZ>);if (pcl::io::loadPCDFile<pcl::PointXYZ> ("second_radius_cloud.pcd", *cloud) == -1){PCL_ERROR("couldn't read file");return 0;}std::cout << "Loaded " << cloud->width * cloud->height<< " data points" << std::endl;// 计算法线pcl::NormalEstimation<pcl::PointXYZ, pcl::Normal> ne;ne.setInputCloud(cloud);pcl::search::KdTree<pcl::PointXYZ>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZ>());ne.setSearchMethod(tree);pcl::PointCloud<pcl::Normal>::Ptr cloud_normals(new pcl::PointCloud<pcl::Normal>);ne.setRadiusSearch(30); // 设置法线估计的半径ne.compute(*cloud_normals);// 法线空间采样pcl::NormalSpaceSampling<pcl::PointXYZ, pcl::Normal> nss;nss.setInputCloud(cloud);nss.setNormals(cloud_normals);nss.setBins(5, 5, 5); // 设置法线空间的bin数量nss.setSample(cloud->size() / 10); // 例如,取原始点云大小的1/10nss.filter(*cloud_downsampled);std::cout << "downsampled cloud size: " << cloud_downsampled->width * cloud_downsampled->height << std::endl;pcl::visualization::PCLVisualizer::Ptr viewer(new pcl::visualization::PCLVisualizer("3D Viewer"));viewer->setBackgroundColor(0, 0, 0); // 设置背景色viewer->addPointCloud<pcl::PointXYZ>(cloud_downsampled, "sample cloud");viewer->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "sample cloud");viewer->initCameraParameters();while (!viewer->wasStopped()){viewer->spinOnce(100);}return 0;
}