RANSAC空间圆拟合实现

news2024/11/18 9:29:20

由初中的几何知识我们可以知道,确定一个三角形至少需要三个不共线的点,因此确定一个三角形的外接圆至少可用三个点。我们不妨假设三个点坐标为P1(x1,y1,z1),P2(x2,y2,z2),P3(x3,y3,z3)。
圆方程的标准形式为:
(xi-x)2+(yi-y)2=R2 (1)
在三维空间坐标系中球的方程为:
(xi-x)2+(yi-y)2+(zi-z)2=R2 (2)
一个空间圆的产生可以看作过该圆心的一个球体,被一个经过该点的平面所截而得到。因此在求取空间圆的时候还应添加平面约束方程,以上述P1、P2、P3三个不共线的点可以确定一个平面,平面方程为:
AX+BY+CZ+D=0 (3)
求解方程(3)时注意技巧,常用两种方式,向量法或者待定系数法。若使用待定系数法则需平面方程的另一种形式,点法式为:
A(x-x0)+B(y-y0)+C(z-z0)=0 (4)
但是笔者建议使用向量的方式,更为简单方便,三点在同一个平面上,以P1为基点,寻找两条向量P12、P13;该平面方程的法向量为n=P12×P13(注:向量的叉乘)。在得到平面的法向量之后,带入其中一个点到方程(4)中即可得到平面约束方程。我们不妨将该约束方程系数设为A1、B1、C1、D1。
在得到约束平面方程之后还须求解空间中圆的方程,在求解圆的方程时也有诸多方式,一种方式为将已知点带入到方程(2)中,利用待定系数法求解,我们仍参以P1为方程的基点,分别将P2,P3带入可以得到如下方程:
(x1-x)2+(y1-y)2+(z1-z)2=R2
(x2-x)2+(y2-y)2+(z2-z)2=R2
(x3-x)2+(y3-y)2+(z3-z)2=R2
整理后可以得到
A2=2(x2-x1) ,B2=2(y2-y1) ,C2=2(z2-z1),D2=-(x12+y12+z12-x22-y22-z22) (5)
A3=2(x3-x1) ,B3=2(y3-y1) ,C3=2(z3-z1),D3=-(x12+y12+z12-x32-y32-z32) (6)
建立系数矩阵:
在这里插入图片描述
求解上述方程即可得到空间圆心坐标(x,y,z)。

RANSAC空间圆拟合的代码实现如下:

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


int main()
{
	pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
	pcl::io::loadPCDFile("circle_3D.pcd", *cloud);

	int iters = 100;					//迭代次数
	float F_thre = 0.1, D_thre = 0.1;	//点到空间圆所在平面的距离阈值,半径误差阈值
	int max_inner_count = 0;			//最大内点个数
	pcl::PointXYZ O;					//圆心
	float R;							//半径
	srand(time(0));						//随机数种子

	for (size_t i = 0; i < iters; i++)	//循环迭代
	{
		int n1 = rand() % cloud->size(), n2 = rand() % cloud->size(), n3 = rand() % cloud->size();
		pcl::PointXYZ p1 = cloud->points[n1], p2 = cloud->points[n2], p3 = cloud->points[n3];	//从点云随机取三个点

		float A1 = (p2.y - p1.y) * (p3.z - p1.z) - (p2.z - p1.z) * (p3.y - p1.y);	//计算空间圆所在平面方程
		float B1 = (p2.z - p1.z) * (p3.x - p1.x) - (p2.x - p1.x) * (p3.z - p1.z);
		float C1 = (p2.x - p1.x) * (p3.y - p1.y) - (p2.y - p1.y) * (p3.x - p1.x);
		float D1 = -(A1 * p1.x + B1 * p1.y + C1 * p1.z);

		float A2 = 2 * (p2.x - p1.x);
		float B2 = 2 * (p2.y - p1.y);
		float C2 = 2 * (p2.z - p1.z);
		float D2 = p1.x * p1.x + p1.y * p1.y + p1.z * p1.z - p2.x * p2.x - p2.y * p2.y - p2.z * p2.z;

		float A3 = 2 * (p3.x - p1.x);
		float B3 = 2 * (p3.y - p1.y);
		float C3 = 2 * (p3.z - p1.z);
		float D3 = p1.x * p1.x + p1.y * p1.y + p1.z * p1.z - p3.x * p3.x - p3.y * p3.y - p3.z * p3.z;

		Eigen::Matrix3f A;
		A << A1, B1, C1, A2, B2, C2, A3, B3, C3;

		Eigen::Vector3f D;
		D << -D1, -D2, -D3;

		auto X = A.inverse() * D;

		pcl::PointXYZ p(X[0], X[1], X[2]);
		float r = (pcl::euclideanDistance(p, p1) + pcl::euclideanDistance(p, p2) + pcl::euclideanDistance(p, p3)) / 3;
		//std::cout << r << std::endl;
		int inner_count = 0;
		for (size_t i = 0; i < cloud->size(); i++)
		{
			float Fi = abs(A1 * cloud->points[i].x + B1 * cloud->points[i].y + C1 * cloud->points[i].z + D1) / sqrt(A1 * A1 + B1 * B1 + C1 * C1);	//点到平面距离
			float Di = abs(pcl::euclideanDistance(p, cloud->points[i]) - r);																		//点到圆心距离于半径差值
			if (Fi < F_thre && Di < D_thre)
				inner_count++;
		}
		if (inner_count > max_inner_count)	//如果本轮迭代内点个数更多,则更新圆心和半径
		{
			max_inner_count = inner_count;
			O = p;
			R = r;
		}
	}
	
	std::cout << O << " " << R << std::endl;

	return 0;
}

参考:在空间三维坐标系下的圆、直线和平面拟合

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/1884199.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

[吃瓜教程]南瓜书第4章决策树

1.决策树的算法原理 从逻辑角度&#xff0c;条件判断语句的组合&#xff1b;从几何角度&#xff0c;根据某种准则划分特征空间&#xff1b; 是一种分治的思想&#xff0c;其最终目的是将样本约分约纯&#xff0c;而划分的核心是在条件的选择或者说是**特征空间的划分标准 ** …

Fooocus模型配置中文教程

很多同学这里不知道该怎么选择。不知道每个模型效果&#xff0c;针对这个整理了一个表格。参考表格就可生成预期效果图。 下载地址&#xff1a; https://download.csdn.net/download/yuanshiren133/89503764

【详解】RV1106移植opencv-mobile库

文章目录 前言一、烧入镜像二、编译项目1.创建项目文件 三、移植四、运行文件五、总结 前言 硬件&#xff1a;瑞芯微Rv1106【Luckfox Pro\Max Pico、网线一根、USB线、串口助手、摄像头 软件&#xff1a;ubuntu 20.4 编译器&#xff1a;arm-rockchip830-linux-uclibcgnueabihf…

Cesium大屏-vue3注册全局组件

1.需求 说明&#xff1a;产品经理要求开发人员在地图大屏上面随意放置组件&#xff0c;并且需要通过数据库更改其组件大小&#xff0c;位置等&#xff1b;适用于大屏组件中场站视角、任意位置标题等。 2.实现 2.1GlobalComponents.vue 说明&#xff1a;containerList可以通…

阿里云物联网应用层开发:第三部分,微信小程序和web客户端实现

文章目录 哔哩哔哩视频教程1、阿里云物联网平台对接微信小程序2、阿里云物联网平台对接web客户端2-1MQTT服务器编写2-2 web端Servlet部分编写 哔哩哔哩视频教程 【阿里云物联网综合开发&#xff0c;STM32ESP8266微信小程序web客户端一篇教程详细讲解】 https://www.bilibili.c…

袋鼠快跳 - 常用网址快捷访问

袋鼠快跳 开源地址&#xff1a;https://github.com/chenbimo/kangaroo-jump 袋鼠快跳&#xff0c;是一个以 简单快捷 为目标的网站快导航 油猴脚本。 本工具的理念就是&#xff0c;用最快的速度访问我们最常用的50个网站。 功能特点 完全免费&#xff0c;以 MIT协议 开源。…

文生图功能介绍

Stable Diffusion WebUI&#xff08;SD WebUI&#xff09;及文生图功能介绍 一、引言 随着人工智能技术的飞速发展&#xff0c;AI绘画作为一种新兴的艺术形式&#xff0c;逐渐走入人们的视野。Stable Diffusion WebUI&#xff08;简称SD WebUI&#xff09;作为AI绘画领域的重…

如何现代的编译和安装内核

前言&#xff1a;本文是在阅读书目时找到了一篇非常高质量的文章。的原文是英文&#xff0c;现在我自己手头翻译了一下&#xff0c;发布到这里。 原文连接&#xff1a;How to compile a Linux kernel in the 21st century | Opensource.com 目录 更新内核的现代方法 安装内…

在线如何快速把图片变小?图片轻松修改大小的3个在线工具

随着现在图片在工作和生活中的广泛使用&#xff0c;在使用图片的时候经常会因为图片太大的问题受到影响&#xff0c;比较简单的一种处理方法可以通过压缩图片的方式来缩小图片大小&#xff0c;那么图片压缩具体该怎么来操作呢&#xff1f;下面就给大家分享几款图片在线压缩工具…

npm安装包报错解决

目录 一&#xff1a;问题回顾 二:问题分析 三&#xff1a;npm降级或者升级 四&#xff1a;npm和node js 关系 一&#xff1a;问题回顾 今天在本地部署一个vue开发的项目&#xff0c;需要在本地看下运行情况&#xff0c;按照常规的操作就是在网站根目录运行npm install 安装…

Android super.img结构及解包和重新组包

Android super.img结构及解包和重新组包 从Android10版本开始&#xff0c;Android系统使用动态分区&#xff0c;system、vendor、 odm等都包含在super.img里面&#xff0c;编译后的最终镜像不再有这些单独的 image&#xff0c;取而代之的是一个总的 super.img. 1. 基础知识 …

LVS-负载均衡

目录 一、概念 二、LVS工作原理 1. ipvs/ipvsadm 2.名词&#xff1a; 三、常用命令 四、工作模式 1.NAT地址转换模式 &#xff08;1&#xff09;工作流程 &#xff08;2&#xff09;特点 &#xff08;3&#xff09;实验过程 a.环境准备&#xff1a; b.修改测试机的…

第二周:计算机网络概述(下)

一、计算机网络性能指标&#xff08;速率、带宽、延迟&#xff09; 1、速率 2、带宽 3、延迟/时延 前面讲分组交换的时候介绍了&#xff0c;有一种延迟叫“传输延迟”&#xff0c;即发送一个报文&#xff0c;从第一个分组的发送&#xff0c;到最后一个分组的发送完成的这段时…

vue3.2及以上 父调子的方法defineExpose定义供父调用的方法及属性

1、定义子类LoginForm&#xff1a; function handleLogin(account, token) {console.log(account,token)}defineExpose({handleLogin,}); 2、父类调用子类组件 const loginFormRef ref(); <LoginForm ref"loginFormRef" />loginFormRef.value.handleLogin(…

仓库管理系统23--用户管理

原创不易&#xff0c;打字不易&#xff0c;截图不易&#xff0c;多多点赞&#xff0c;送人玫瑰&#xff0c;留有余香&#xff0c;财务自由明日实现 1、创建用户管理的用户控件 <UserControl x:Class"West.StoreMgr.View.UserInfoView"xmlns"http://schemas.…

SSE代替轮询?

什么是 SSE SSE&#xff08;Server-Sent Events&#xff0c;服务器发送事件&#xff09;&#xff0c;为特定目的而扩展的 HTTP 协议&#xff0c;用于实现服务器向客户端推送实时数据的单向通信。如果连接断开&#xff0c;浏览器会自动重连&#xff0c;传输的数据基于文本格式。…

C语言入门-指针和数组5

指针和地址 地址 地址是内存中一个特定位置的标识符。每个内存位置都有一个唯一的地址&#xff0c;用于存储数据。这些地址通常表示为十六进制数。 物理地址&#xff1a;硬件层次上的实际内存地址。逻辑地址&#xff1a;程序运行时使用的地址&#xff0c;由操作系统管理。 …

Edge浏览器添加新标签页网址为 百度 搜索

默认不能直接设置&#xff0c;需要安装New Tab Change插件 安装拓展插件 url 这里点击获取即可&#xff08;我已经安装过了&#xff09; 某些扩展会更改浏览器设置&#xff0c;例如默认搜索引擎、新选项卡页和其他类型的网站数据。 为了防止扩展更改在安装时设置的首选项Micr…

MQ:RabbitMQ

同步和异步通讯 同步通讯: 需要实时响应,时效性强 耦合度高 每次增加功能都要修改两边的代码 性能下降 需要等待服务提供者的响应,如果调用链过长则每次响应时间需要等待所有调用完成 资源浪费 调用链中的每个服务在等待响应过程中,不能释放请求占用的资源,高并发场景下…

来聊聊Redis定期删除策略的设计与实现

写在文章开头 我们都知道redis通过主线程完成内存数据库的指令操作&#xff0c;由于只有一个线程负责核心业务流程&#xff0c;所以对于每一个操作都要求尽可能达到尽可能的高效迅速&#xff0c;而本文就基于源码来聊聊redis的定期删除策略的设计与实现。 Hi&#xff0c;我是 …