点云边界提取方法总结

2023-11-10

原始点云:
在这里插入图片描述

经纬线扫描法

经纬线扫描法的基本思想是:将经过坐标变换后的二维数据点集按照 x值的大小排序后存储在一个链表中,在二维平面建立点集的最小包围盒并分别计算出 x 和 y 的最大、最小值;令经线从 x 的最小值开始,取步长为dx,在 x 的取值范围内分别计算出每根经线的最大和最小 y 值,并将它们的索引值放在一个新建的链表中,扫描完整个 x 区间;同样的方法扫描纬线,最后将重复的索引值删掉。
该算法的运算速度受点云点数以及搜索分辨率影响。下面代码在点云形状近似为平面,且在xoy平面上投影面积最大时表现较好,其他情况需要修改。

#include <iostream>
#include <algorithm>
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>

void BoundaryExtraction(pcl::PointCloud<pcl::PointXYZ>::Ptr cloud, pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_boundary, int resolution)
{
	pcl::PointXYZ px_min = *std::min_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.x < pt2.x; });
	pcl::PointXYZ px_max = *std::max_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.x < pt2.x; });

	float delta_x = (px_max.x - px_min.x) / resolution;
	float min_y = INT_MAX, max_y = -INT_MAX;
	std::vector<int> indexs_x(2 * resolution + 2);
	std::vector<std::pair<float, float>> minmax_x(resolution + 1, { INT_MAX,-INT_MAX });
	for (size_t i = 0; i < cloud->size(); ++i)
	{
		int id = (cloud->points[i].x - px_min.x) / delta_x;
		if (cloud->points[i].y < minmax_x[id].first)
		{
			minmax_x[id].first = cloud->points[i].y;
			indexs_x[id] = i;
		}
		else if (cloud->points[i].y > minmax_x[id].second)
		{
			minmax_x[id].second = cloud->points[i].y;
			indexs_x[id + resolution + 1] = i;
		}
	}

	pcl::PointXYZ py_min = *std::min_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.y < pt2.y; });
	pcl::PointXYZ py_max = *std::max_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.y < pt2.y; });

	float delta_y = (py_max.y - py_min.y) / resolution;
	float min_x = INT_MAX, max_x = -INT_MAX;
	std::vector<int> indexs_y(2 * resolution + 2);
	std::vector<std::pair<float, float>> minmax_y(resolution + 1, { INT_MAX,-INT_MAX });
	for (size_t i = 0; i < cloud->size(); ++i)
	{
		int id = (cloud->points[i].y - py_min.y) / delta_y;
		if (cloud->points[i].x < minmax_y[id].first)
		{
			minmax_y[id].first = cloud->points[i].x;
			indexs_y[id] = i;
		}
		else if (cloud->points[i].x > minmax_y[id].second)
		{
			minmax_y[id].second = cloud->points[i].x;
			indexs_y[id + resolution + 1] = i;
		}
	}

	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_xboundary(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::copyPointCloud(*cloud, indexs_x, *cloud_xboundary);
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_yboundary(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::copyPointCloud(*cloud, indexs_y, *cloud_yboundary);
	*cloud_boundary = *cloud_xboundary + *cloud_yboundary;
}

int main(int argc, char* argv[])
{
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>),;
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_boundary(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("test.pcd", *cloud);
	
	BoundaryExtraction(cloud, cloud_boundary, 200);
	pcl::io::savePCDFile("boundary.pcd", *cloud_boundary);

	return EXIT_SUCCESS;
}

提取结果:
在这里插入图片描述

网格划分法

网格划分法分为三个步骤:(1)网格划分(2)寻找边界网格(3)提取边界线。首先建立数据点集的最小包围盒,并用给定间隔的矩形网格将其分割;然后找出边界网格,将这些边界网格按顺序连接起来形成一条由边界网格组成的“粗边界”;再对每个边界网格按照某种规则判断其内的点是否为边界点,连接初始边界,并对此边界线进一步处理使其光滑。
在这里插入图片描述
在这里插入图片描述
该算法的运算速度受点云点数以及搜索分辨率影响。下面代码在点云形状近似为平面,且在xoy平面上投影面积最大时表现较好,其他情况需要修改。

#include <iostream>
#include <vector>
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/common/transforms.h>

void BoundaryExtraction( pcl::PointCloud<pcl::PointXYZ>::Ptr cloud,  pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_boundary, int resolution)
{
	pcl::PointXYZ px_min = *std::min_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.x < pt2.x; });
	pcl::PointXYZ px_max = *std::max_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.x < pt2.x; });
	pcl::PointXYZ py_min = *std::min_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.y < pt2.y; });
	pcl::PointXYZ py_max = *std::max_element(cloud->begin(), cloud->end(), [](pcl::PointXYZ pt1, pcl::PointXYZ pt2) {return pt1.y < pt2.y; });
	float x_min = px_min.x, x_max = px_max.x, y_min = py_min.y, y_max = py_max.y;

	float L = sqrt((x_max - x_min)*(y_max - y_min) / resolution);
	int x_num = (x_max - x_min) / L + 1;
	int y_num = (y_max - y_min) / L + 1;

	std::vector<std::vector< pcl::PointCloud<pcl::PointXYZ>::Ptr>> uv(x_num + 1, std::vector< pcl::PointCloud<pcl::PointXYZ>::Ptr>(y_num + 1));
	for (int i = 0; i <= x_num; ++i)
	{
		for (int j = 0; j <= y_num; ++j)
		{
			 pcl::PointCloud<pcl::PointXYZ>::Ptr ptcloud(new  pcl::PointCloud<pcl::PointXYZ>);
			 uv[i][j] = ptcloud;
		}
	}
	for (int i = 0; i < cloud->size(); ++i)
	{
		int x_id = (cloud->points[i].x - x_min) / L;
		int y_id = (cloud->points[i].y - y_min) / L;
		uv[x_id][y_id]->push_back(cloud->points[i]);
	}

	std::vector<std::vector<bool>> boundary_index(x_num + 1, std::vector<bool>(y_num + 1, false));
	for (int i = 0; i <= x_num; ++i)
	{
		if (uv[i][0]->size())
			boundary_index[i][0] = true;
		if (uv[i][y_num]->size())
			boundary_index[i][y_num] = true;
	}		
	for (int i = 0; i <= y_num; ++i)
	{
		if(uv[0][i]->size())
			boundary_index[0][i] = true;
		if (uv[x_num][i]->size())
			boundary_index[x_num][i] = true;
	}
	for (int i = 1; i < x_num - 1; ++i)
	{
		for (int j = 1; j < y_num - 1; ++j)
		{
			if (uv[i][j]->size())
			{
				if (!uv[i - 1][j - 1]->size() || !uv[i][j - 1]->size() || !uv[i + 1][j - 1]->size() || !uv[i][j - 1]->size() ||
					!uv[i + 1][j - 1]->size() || !uv[i + 1][j]->size() || !uv[i + 1][j + 1]->size() || !uv[i][j + 1]->size())
					boundary_index[i][j] = true;
			}
		}
	}

	for (int i = 0; i <= x_num; ++i)
	{
		for (int j = 0; j <= y_num; ++j)
		{
			if (boundary_index[i][j])
			{
				 pcl::PointCloud<pcl::PointXYZ>::Ptr ptcloud(new  pcl::PointCloud<pcl::PointXYZ>);
				ptcloud = uv[i][j];
				Eigen::Vector4f cloud_centroid;
				pcl::compute3DCentroid(*ptcloud, cloud_centroid);
				cloud_boundary->push_back(pcl::PointXYZ(cloud_centroid(0), cloud_centroid(1), cloud_centroid(2)));
			}
		}
	}
}

int main(int argc, char* argv[])
{
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new  pcl::PointCloud<pcl::PointXYZ>);
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_boundary(new  pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("test.pcd", *cloud);
	
	BoundaryExtraction(cloud, cloud_boundary, 200*200);
	pcl::io::savePCDFile("boundary.pcd", *cloud_boundary);

	return EXIT_SUCCESS;
}

提取结果:
在这里插入图片描述

法线估计法

可参考PCL:点云数据基于法线的边界提取(从最初的法线估计理论推导到最终的边界提取)
该算法的运算速度受点云点数以及搜索半径影响。

#include <iostream>
#include <pcl/io/pcd_io.h>
#include <pcl/features/normal_3d.h>
#include <pcl/features/boundary.h>


int main(int argc, char* argv[])
{

	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("test.pcd", *cloud);

	pcl::PointCloud<pcl::Normal>::Ptr normals(new pcl::PointCloud<pcl::Normal>);
	pcl::PointCloud<pcl::Boundary> boundaries;
	pcl::BoundaryEstimation<pcl::PointXYZ, pcl::Normal, pcl::Boundary> est;
	pcl::search::KdTree<pcl::PointXYZ>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZ>());

	pcl::NormalEstimation<pcl::PointXYZ, pcl::Normal> normEst;
	normEst.setInputCloud(cloud);
	normEst.setSearchMethod(tree);
	normEst.setRadiusSearch(3);  //法向估计的半径  3
	//normEst.setKSearch(9);  //法向估计的点数
	normEst.compute(*normals);

	est.setInputCloud(cloud);
	est.setInputNormals(normals);
	est.setSearchMethod(tree);
	est.setKSearch(500);  //一般这里的数值越高,最终边界识别的精度越好
	est.compute(boundaries);

	pcl::PointCloud<pcl::PointXYZ>::Ptr boundPoints(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::PointCloud<pcl::PointXYZ> noBoundPoints;
	int countBoundaries = 0;
	for (int i = 0; i < cloud->size(); i++)
	{
		uint8_t x = (boundaries.points[i].boundary_point);
		int a = static_cast<int>(x); //该函数的功能是强制类型转换
		if (a == 1)
		{
			(*boundPoints).push_back(cloud->points[i]);
			countBoundaries++;
		}
		else
			noBoundPoints.push_back(cloud->points[i]);
	}

	pcl::io::savePCDFile("boundary.pcd", *boundPoints);

	return EXIT_SUCCESS;
}

提取结果:
在这里插入图片描述

alpha shapes算法

原理参考:alpha shapes提取平面点云边界点
该算法的运算速度受点云点数以及搜索半径影响。
手写实现:在alpha-shape计算二维点云边界–c++实现基础上精简优化

#include <iostream>
#include <pcl/io/pcd_io.h>
#include <pcl/common/distances.h>
#include <pcl/kdtree/kdtree_flann.h>


int main()
{
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("test.pcd", *cloud);
	for (auto& point : *cloud)
		point.z = 0.0;

	pcl::PointCloud<pcl::PointXYZ>::Ptr boundary(new pcl::PointCloud<pcl::PointXYZ>);
	float alpha = 10.0;
	pcl::KdTreeFLANN<pcl::PointXYZ>kdtree;
	kdtree.setInputCloud(cloud);
	std::vector<int> indices;
	std::vector<float>dist;
	for (size_t i = 0; i < cloud->size(); ++i)
	{
		pcl::PointXYZ p = cloud->points[i];
		kdtree.radiusSearch(cloud->points[i], 2*alpha, indices, dist); //搜索在点云中到p点距离小于2*alpha的点
		for (size_t j = 1; j < indices.size(); ++j)
		{
			pcl::PointXYZ p1 = cloud->points[indices[j]];
			float s_2 = std::pow((p.x - p1.x), 2) + std::pow((p.y - p1.y), 2);
			float h = std::pow((alpha * alpha / s_2 - 0.25), 0.5);
			float x2 = p.x + 0.5 * (p1.x - p.x) - h * (p1.y - p.y);
			float y2 = p.y + 0.5 * (p1.y - p.y) - h * (p.x - p1.x);
			float x3 = p.x + 0.5 * (p1.x - p.x) + h * (p1.y - p.y);
			float y3 = p.y + 0.5 * (p1.y - p.y) + h * (p.x - p1.x);
			pcl::PointXYZ p2(x2, y2, 0.0);
			pcl::PointXYZ p3(x3, y3, 0.0);

			//计算邻域内除了p1之外的点到p2、p3的距离
			std::vector<bool> distp2_bool, distp3_bool;
			int count = 0;
			for (size_t k = 1; k < indices.size(); ++k)
			{
				pcl::PointXYZ p_other = cloud->points[indices[k]];
				if (k != j)
				{
					++count;
					//比较剩余的点到p2、p3距离与r的大小
					if (pcl::euclideanDistance(p_other, p2) > alpha)
						distp2_bool.push_back(true);
					if (pcl::euclideanDistance(p_other, p3) > alpha)
						distp3_bool.push_back(true);
				}
			}
			//如果邻域内所有点到p2,或者p3的距离均大于alpha,则表明p是边界点,退出循环
			if (count == distp2_bool.size() || count == distp3_bool.size())
			{
				boundary->push_back(cloud->points[i]);
				break;
			}
		}
	}

	pcl::io::savePCDFile("boundary.pcd", *boundary);

	return EXIT_SUCCESS;
}

调用api:

#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/surface/concave_hull.h>


int main(int argc, char* argv[])
{
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("test.pcd", *cloud);

	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_hull(new pcl::PointCloud<pcl::PointXYZ>); 
	pcl::ConcaveHull<pcl::PointXYZ> chull;
	chull.setInputCloud(cloud); //输入点云为投影后的点云
	chull.setAlpha(10); //设置aLpha值
	chull.reconstruct(*cloud_hull);
	pcl::io::savePCDFile("boundary.pcd", *cloud_hull);

	return EXIT_SUCCESS;
}

提取结果:
在这里插入图片描述

参考文献:[1]杨雪娇. 点云的边界提取及角点检测算法研究[D].哈尔滨工程大学,2010.

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

点云边界提取方法总结 的相关文章

  • Dynamics 365详解

    什么是Dynamics 365 Dynamics 365是微软公司推出的一款企业资源计划 ERP 和客户关系管理 CRM 软件 它是微软旗下的云计算平台Azure上的一项服务 能够在多个设备和平台上运行 Dynamics 365结合了ERP
  • leetcode分类刷题:队列(Queue)(三、优先队列用于归并排序)

    1 当TopK问题出现在多个有序序列中时 就要用到归并排序的思想了 2 将优先队列初始化为添加多个有序序列的首元素的形式 再循环K次优先队列的出队和出队元素对应序列下个元素的入队 就能得到TopK的元素了 3 这些题目好像没有TopK 大用

随机推荐

  • Python--Email

    邮件定义 电子邮件 消息由头域 统称消息头 以及后面可选的消息体组成 根据 RFC 2822 唯一需要的消息标题只有发送日期字段和发送地址字段 即 Date 和 From MAIL FROM RCPT TO DATA 1 电子邮件系统组件和
  • Jbox2D入门学习一物理世界及最简单的物体创建

    这周末无聊 翻到舍友有本游戏开发的书 就浏览了遍 因为对游戏以前其实没接触过 可能就简单知道一边游戏绘图和逻辑以及游戏框架绘制的简单概念 这次主要是看了下游戏开发中最常用和基础的一个物理引擎 Box2D 对于这个引擎 可能说最好表现和表达出
  • Spring Boot彩色日志配置

    在spring boot中使用彩色日志 spring boot是默认支持彩色日志的 但是由于我又添加了自己的logback日志配置文件 然后就没有了彩色日志 经过一番搜索大法找到了一个完美还原spring boot的彩色日志 下面是logb
  • 华为OD机试 - 数字加减游戏(Java)

    题目描述 小明在玩一个数字加减游戏 只使用加法或者减法 将一个数字s变成数字t 每个回合 小明可以用当前的数字加上或减去一个数字 现在有两种数字可以用来加减 分别为a b a b 其中b没有使用次数限制 请问小明最少可以用多少次a 才能将数
  • 技术流

    如果让你向别人推荐十部电影 你会推荐哪十部 这是在知乎上被浏览过 2000 万的问题 一共有 5036 个回答 花 5 分钟读完这篇文章 可以帮你节省99 的找电影时间 全是知乎上最值得推荐的电影 并最后让你获取到几百万人收藏的超级电影名单
  • Spring 之 jwt,过滤器,拦截器,aop,监听器,参数校验

    Spring 之 jwt 过滤器 拦截器 aop 监听器 一 jwt编写 1 1 pom 1 2 JwtUtils 1 3 注意 1 4 用法 1 5 jwt实战 1 5 1 过滤器判断jwt 1 5 2 过滤器注入Bean 1 5 3 t
  • under-approximation & over-approximation

    Under approximation and over approximation are concepts often used in the context of formal methods a field that applies
  • IDM403解决办法之一:和百度网盘和平共处

    不知道有没有小伙伴和我一样 对百度云恨之入骨 在千兆时代 它的下载速度竟然还能再100KB左右 于是 我在网络上找了各种方法 最好用的是用IDM下载 不过有的时候 IDM也会失灵 特别是针对百度网盘里的大文件 以前300MB以上 这次我下载
  • git代码使用空格缩进

    1 idea设置缩进符为空格 Java 代码 golang 代码 2 设置提交仓库时的空格处理 否则 golang 代码为了减少文件大小 可能会把空格缩进改为制表符 设置当前仓库配置 git config core whitespace t
  • sublime text3 智能提示

    1 启动编辑器 编写代码方法 发现只有html方面的提示 并没有函数方面的提示 2 点击 preferences package control 准备插件包各项操作 3 点击下拉框中的 install package 安装插件包 4 在输入
  • 如何提高一个研发团队的“代码速度”?

    阿里妹导读 Code Velocity 代码速度 体现了一个研发团队快速响应业务需求的能力 如果做得好 代码从commit到上线可能平均只需要两三天时间 甚至连紧急发布都不怎么需要了 今天 蚂蚁金服国际事业群技术风险部研究员南门 将和大家聊
  • linux的oops界面,Linux编程时遇到Oops提示该如何排查?

    用于表示函数的调用关系 通过这段信息我们可以知道 函数的整个执行流程 知道它的函数调用关系 最后整理出来的函数执行流程如下 本文引用地址 http www eepw com cn article 201811 395128 htm 从中我们
  • 基于LiDAR的目标检测算法

    自动驾驶中 激光雷达点云如何做特征表达 基于激光雷达点云 lidar 的目标检测方法之BEV 基于激光雷达点云 lidar 的目标检测方法之camera range view 基于激光雷达点云 lidar 的目标检测方法之point wis
  • VueTouchKeyboard——一个模拟键盘

    功能需求 封装一个带有设计好的样式的输入组件 输入方式为模拟的数字键盘 键盘组件为VueTouchKeyboard 下载方式如下 npm install vue touch keyboard save image png 封装的输入组件 模
  • 素数的求解方法:

    一 朴素判断素数算法 就判断素数而言 事实上是非常简单的了 根据定义 判断一个整数n是否是素数 只需要去判断在整数区间 2 n 1 之内 是否具有某个数m 使得n m 0 代码可以这么写 int isPrime int n int i fo
  • Arduino IDE 烧录 ESP8266教程

    Arduino IDE for ESP8266教程 原出处 http www windworkshop cn p 758 ESP8266是现在性价比不错的Wifi模块 用了一块ESP8266 01之后感觉还行 用在数据采集器上表现还是不错的
  • IDEA中build path导入外包jar 包的方法

    小编今天需要在IDEA中导入外部jar包 由于在eclipse中可以在外部jar包上直接右键 build path 然后点击add to Build Path 就成功导入了 可惜IDEA并不是如此 有点小苦恼 接下来 小编就带大家操作 ID
  • 分解命令行字符串为argc和argv

    有时候需要用空格把一个命令行参数字符串分解为参数个数和参数指针 就是常见的c语言main 函数入口argc argv 这里采用strtok 函数可以很方便的做到 char strtok char str const char delim 用
  • 关于取消开机bois密码的方法

    有时在不小心设置了bois密码 导致每次开机都需要输入密码 显得十分麻烦 如何取消bois密码呢 经过查阅各种资料 有两种2情况 分别是 知道bois密码 和不知道bois密码 大部分的文章中都讲到的是忘记bois密码的方式 这边将一下 知
  • 点云边界提取方法总结

    目录 经纬线扫描法 网格划分法 法线估计法 alpha shapes算法 原始点云 经纬线扫描法 经纬线扫描法的基本思想是 将经过坐标变换后的二维数据点集按照 x值的大小排序后存储在一个链表中 在二维平面建立点集的最小包围盒并分别计算出 x