C++ 计算凸包点的最小旋转矩形

RotateRect.h


#include <vector>

/**
* @brief 计算点集最小旋转外接矩形
*/
class RotateRect {
public:
	enum { CALIPERS_MAXHEIGHT = 0, CALIPERS_MINAREARECT = 1, CALIPERS_MAXDIST = 2 };
	struct Point {
		float x, y;
	};
	using Points = std::vector<Point>;
	struct Sizef {
		float width, height;
	};
	struct Rect {
		float left, top, right, bottom;
	};
	explicit RotateRect(const Points& pts);

	const Point& center() const;
	const Sizef& size() const;
	float angle() const;
public:
	void update();
	Points toPoints() const;
	Rect getOutLine() const;

private:
	double crossProduct(const Point& A, const Point& B, const Point& C) const;
	double calculateDistance(const Point& A, const Point& B) const;
	Points convexHull() const;
	/* we will use usual cartesian coordinates */
	void rotatingCalipers(const Point* pts, int n, int mode, float* out) const;

private:
	const Points& m_inputs;
	Point m_center;
	//! returns width and height of the rectangle
	Sizef m_size;
	//! returns the rotation angle. When the angle is 0, 90, 180, 270 etc., the rectangle becomes an up-right rectangle.
	float m_angle;
};

 RotateRect.cpp


#include "RotateRect.h"
#include <cmath>
#include <assert.h>

#ifndef CV_PI
#define CV_PI   3.1415926535897932384626433832795
#endif

RotateRect::RotateRect(const RotateRect::Points& pts) : m_inputs(pts), m_angle(0), m_size{ 0, 0 }, m_center{ 0, 0 } {}

/*
 * @brief 计算BA与CA之间面积
 * @param [in] A 点
 * @param [in] B 点
 * @param [in] C 点
 * @return
	BA x CA > 0, C在B的逆时针方向,B在C右边
	BA x CA < 0, C在B的顺时针方向,B在C左边
	BA x CA = 0, C和B共线,可能正方向,也可能反方向
 */
double RotateRect::crossProduct(const RotateRect::Point& A, const RotateRect::Point& B, const RotateRect::Point& C) const {
	return (B.x - A.x) * (C.y - A.y) - (C.x - A.x) * (B.y - A.y);
}


/*
 * @brief 找到x方向最小点为起始点,查找所有点最右边点,到起始点结束
 * @return 所有凸包点,逆时针方向
 */
RotateRect::Points RotateRect::convexHull() const {
	const auto& points = m_inputs;
	int n = points.size();
	if (n <= 3) {
		return points;
	}

	std::vector<Point> hull;
	int l = 0;
	for (int i = 1; i < n; i++) {
		if (points[i].x < points[l].x) {
			l = i;
		}
	}

	int p = l, q;
	do {
		hull.push_back(points[p]);
		q = (p + 1) % n;
		for (int i = 0; i < n; i++) {
			if (crossProduct(points[p], points[i], points[q]) > 0) {
				q = i;
			}
		}
		p = q;
	} while (p != l);

	return hull;
}

double RotateRect::calculateDistance(const RotateRect::Point& A, const RotateRect::Point& B) const {
	return std::sqrt(std::pow(B.x - A.x, 2) + std::pow(B.y - A.y, 2));
}

void rotate90CCW(const RotateRect::Point& in, RotateRect::Point& out)
{
	out.x = -in.y;
	out.y = in.x;
}

void rotate90CW(const RotateRect::Point& in, RotateRect::Point& out)
{
	out.x = in.y;
	out.y = -in.x;
}

void rotate180(const RotateRect::Point& in, RotateRect::Point& out)
{
	out.x = -in.x;
	out.y = -in.y;
}

/* return true if first vector is to the right (clockwise) of the second */
bool firstVecIsRight(const RotateRect::Point& vec1, const RotateRect::Point& vec2)
{
	RotateRect::Point tmp;
	rotate90CW(vec1, tmp);
	return tmp.x * vec2.x + tmp.y * vec2.y < 0;
}

/*
 * @brief 计算凸包点的最小旋转矩形
 */
 /* we will use usual cartesian coordinates */
void RotateRect::rotatingCalipers(const RotateRect::Point* points, int n, int mode, float* out) const
{
	using Point = RotateRect::Point;
	float minarea = FLT_MAX;
	float max_dist = 0;
	char buffer[32] = {};
	int i, k;
	std::vector<float> abuf(n * 3);
	float* inv_vect_length = abuf.data();
	Point* vect = (Point*)(inv_vect_length + n);
	int left = 0, bottom = 0, right = 0, top = 0;
	int seq[4] = { -1, -1, -1, -1 };
	Point rot_vect[4];

	/* rotating calipers sides will always have coordinates
	 (a,b) (-b,a) (-a,-b) (b, -a)
	 */
	 /* this is a first base vector (a,b) initialized by (1,0) */
	float orientation = 0;
	float base_a;
	float base_b = 0;

	float left_x, right_x, top_y, bottom_y;
	Point pt0 = points[0];

	left_x = right_x = pt0.x;
	top_y = bottom_y = pt0.y;

	for (i = 0; i < n; i++)
	{
		double dx, dy;

		if (pt0.x < left_x)
			left_x = pt0.x, left = i;

		if (pt0.x > right_x)
			right_x = pt0.x, right = i;

		if (pt0.y > top_y)
			top_y = pt0.y, top = i;

		if (pt0.y < bottom_y)
			bottom_y = pt0.y, bottom = i;

		Point pt = points[(i + 1) & (i + 1 < n ? -1 : 0)];

		dx = pt.x - pt0.x;
		dy = pt.y - pt0.y;

		vect[i].x = (float)dx;
		vect[i].y = (float)dy;
		inv_vect_length[i] = (float)(1. / std::sqrt(dx * dx + dy * dy));

		pt0 = pt;
	}

	// find convex hull orientation
	{
		double ax = vect[n - 1].x;
		double ay = vect[n - 1].y;

		for (i = 0; i < n; i++)
		{
			double bx = vect[i].x;
			double by = vect[i].y;

			double convexity = ax * by - ay * bx;

			if (convexity != 0)
			{
				orientation = (convexity > 0) ? 1.f : (-1.f);
				break;
			}
			ax = bx;
			ay = by;
		}
		assert(orientation != 0);
	}
	base_a = orientation;

	/*****************************************************************************************/
	/*                         init calipers position                                        */
	seq[0] = bottom;
	seq[1] = right;
	seq[2] = top;
	seq[3] = left;
	/*****************************************************************************************/
	/*                         Main loop - evaluate angles and rotate calipers               */

	/* all of edges will be checked while rotating calipers by 90 degrees */
	for (k = 0; k < n; k++)
	{
		/* number of calipers edges, that has minimal angle with edge */
		int main_element = 0;

		/* choose minimum angle between calipers side and polygon edge by dot product sign */
		rot_vect[0] = vect[seq[0]];
		rotate90CW(vect[seq[1]], rot_vect[1]);
		rotate180(vect[seq[2]], rot_vect[2]);
		rotate90CCW(vect[seq[3]], rot_vect[3]);
		for (i = 1; i < 4; i++)
		{
			if (firstVecIsRight(rot_vect[i], rot_vect[main_element]))
				main_element = i;
		}

		/*rotate calipers*/
		{
			//get next base
			int pindex = seq[main_element];
			float lead_x = vect[pindex].x * inv_vect_length[pindex];
			float lead_y = vect[pindex].y * inv_vect_length[pindex];
			switch (main_element)
			{
			case 0:
				base_a = lead_x;
				base_b = lead_y;
				break;
			case 1:
				base_a = lead_y;
				base_b = -lead_x;
				break;
			case 2:
				base_a = -lead_x;
				base_b = -lead_y;
				break;
			case 3:
				base_a = -lead_y;
				base_b = lead_x;
				break;
			default:
				assert("main_element should be 0, 1, 2 or 3" && false);
			}
		}
		/* change base point of main edge */
		seq[main_element] += 1;
		seq[main_element] = (seq[main_element] == n) ? 0 : seq[main_element];

		switch (mode)
		{
		case RotateRect::CALIPERS_MAXHEIGHT:
		{
			/* now main element lies on edge aligned to calipers side */

			/* find opposite element i.e. transform  */
			/* 0->2, 1->3, 2->0, 3->1                */
			int opposite_el = main_element ^ 2;

			float dx = points[seq[opposite_el]].x - points[seq[main_element]].x;
			float dy = points[seq[opposite_el]].y - points[seq[main_element]].y;
			float dist;

			if (main_element & 1)
				dist = (float)fabs(dx * base_a + dy * base_b);
			else
				dist = (float)fabs(dx * (-base_b) + dy * base_a);

			if (dist > max_dist)
				max_dist = dist;
		}
		break;
		case RotateRect::CALIPERS_MINAREARECT:
			/* find area of rectangle */
		{
			float height;
			float area;

			/* find vector left-right */
			float dx = points[seq[1]].x - points[seq[3]].x;
			float dy = points[seq[1]].y - points[seq[3]].y;

			/* dotproduct */
			float width = dx * base_a + dy * base_b;

			/* find vector left-right */
			dx = points[seq[2]].x - points[seq[0]].x;
			dy = points[seq[2]].y - points[seq[0]].y;

			/* dotproduct */
			height = -dx * base_b + dy * base_a;

			area = width * height;
			if (area <= minarea)
			{
				float* buf = (float*)buffer;

				minarea = area;
				/* leftist point */
				((int*)buf)[0] = seq[3];
				buf[1] = base_a;
				buf[2] = width;
				buf[3] = base_b;
				buf[4] = height;
				/* bottom point */
				((int*)buf)[5] = seq[0];
				buf[6] = area;
			}
		}
		break;
		}                       /*switch */
	}                           /* for */

	switch (mode)
	{
	case RotateRect::CALIPERS_MINAREARECT:
	{
		float* buf = (float*)buffer;

		float A1 = buf[1];
		float B1 = buf[3];

		float A2 = -buf[3];
		float B2 = buf[1];

		float C1 = A1 * points[((int*)buf)[0]].x + points[((int*)buf)[0]].y * B1;
		float C2 = A2 * points[((int*)buf)[5]].x + points[((int*)buf)[5]].y * B2;

		float idet = 1.f / (A1 * B2 - A2 * B1);

		float px = (C1 * B2 - C2 * B1) * idet;
		float py = (A1 * C2 - A2 * C1) * idet;

		out[0] = px;
		out[1] = py;

		out[2] = A1 * buf[2];
		out[3] = B1 * buf[2];

		out[4] = A2 * buf[4];
		out[5] = B2 * buf[4];
	}
	break;
	case RotateRect::CALIPERS_MAXHEIGHT:
	{
		out[0] = max_dist;
	}
	break;
	}
}

void RotateRect::update()
{
	using Point = RotateRect::Point;
	std::vector<Point> hull;
	Point out[3];
	hull = convexHull();

	int n = hull.size();
	const Point* hpoints = &hull[0];

	if (n > 2)
	{
		rotatingCalipers(hpoints, n, RotateRect::CALIPERS_MINAREARECT, (float*)out);
		m_center.x = out[0].x + (out[1].x + out[2].x) * 0.5f;
		m_center.y = out[0].y + (out[1].y + out[2].y) * 0.5f;
		m_size.width = (float)std::sqrt((double)out[1].x * out[1].x + (double)out[1].y * out[1].y);
		m_size.height = (float)std::sqrt((double)out[2].x * out[2].x + (double)out[2].y * out[2].y);
		m_angle = (float)atan2((double)out[1].y, (double)out[1].x);
	}
	else if (n == 2)
	{
		m_center.x = (hpoints[0].x + hpoints[1].x) * 0.5f;
		m_center.y = (hpoints[0].y + hpoints[1].y) * 0.5f;
		double dx = hpoints[1].x - hpoints[0].x;
		double dy = hpoints[1].y - hpoints[0].y;
		m_size.width = (float)std::sqrt(dx * dx + dy * dy);
		m_size.height = 0;
		m_angle = (float)atan2(dy, dx);
	}
	else
	{
		if (n == 1)
			m_center = hpoints[0];
	}

	m_angle = (float)(m_angle * 180 / CV_PI);
}


RotateRect::Points RotateRect::toPoints() const
{
	RotateRect::Points pt(4);

	double _angle = m_angle * CV_PI / 180.;
	float b = (float)cos(_angle) * 0.5f;
	float a = (float)sin(_angle) * 0.5f;

	pt[0].x = m_center.x - a * m_size.height - b * m_size.width;
	pt[0].y = m_center.y + b * m_size.height - a * m_size.width;
	pt[1].x = m_center.x + a * m_size.height - b * m_size.width;
	pt[1].y = m_center.y - b * m_size.height - a * m_size.width;
	pt[2].x = 2 * m_center.x - pt[0].x;
	pt[2].y = 2 * m_center.y - pt[0].y;
	pt[3].x = 2 * m_center.x - pt[1].x;
	pt[3].y = 2 * m_center.y - pt[1].y;
	return pt;
}

const RotateRect::Point& RotateRect::center() const {
	return m_center;
}

const RotateRect::Sizef& RotateRect::size() const {
	return m_size;
}

float RotateRect::angle() const {
	return m_angle;
}


RotateRect::Rect RotateRect::getOutLine() const {
	if (m_inputs.empty())
		return { 0, 0, 0, 0 };
	using Number = std::numeric_limits<float>;
	RotateRect::Rect rect{ Number::max(), Number::max(), Number::min(), Number::min() };
	for (const auto& point : m_inputs) {
		if (point.x < rect.left)
			rect.left = point.x;
		if (point.x > rect.right)
			rect.right = point.x;

		if (point.y < rect.top)
			rect.top = point.y;
		if (point.y > rect.bottom)
			rect.bottom = point.y;
	}
	return rect;
}

main.cpp


RotateRect::Points myPoints = { {233, 86}, {322, 106}, {214, 154}, {307, 176}, {286, 209}, {331, 183}, {346, 319}, {392, 294}, {356, 346}, {419, 311}, {419, 311}, {778, 1031}, {840, 995} };
    RotateRect myRotateRect(myPoints);
    myRotateRect.update();

 


创作不易,小小的支持一下吧!

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

如若内容造成侵权/违法违规/事实不符,请联系我们进行投诉反馈qq邮箱809451989@qq.com,一经查实,立即删除!

相关文章

FlinkCDC 3.1.0 与 Flink 1.18.0 安装及使用 Mysql To Doris 整库同步,使用 pipepline连接器

cd flink-cdc-3.1.0 bin/flink-cdc.sh 会用到 linux的系统环境变量&#xff08;vim /etc/profile配置&#xff09;&#xff0c;使用环境变量 FLINK_HOME flinkcdc & flink 安装及使用&#xff1a; 1、flink-cdc-3.1.0/lib/ 内容如下&#xff1a; 2、flink-cdc-3.1.0/mysql…

硅谷裸机云服务器:定义、特点与应用

硅谷&#xff0c;作为全球科技创新的重镇&#xff0c;一直是科技领域的风向标。近年来&#xff0c;随着云计算技术的飞速发展&#xff0c;硅谷的科技公司们也在不断探索和创新&#xff0c;以满足日益增长的计算需求。其中&#xff0c;裸机云服务器作为一种新兴的计算资源交付方…

新手下白对Latex下手啦!

第一次使用latex&#xff0c;浅浅地记录一下子吧。 首先我们一般会下载一个latex模板&#xff0c;如果想知道咋下载&#xff0c;评论去告诉俺哟&#xff01; 新手小白首先要看懂结构&#xff0c;不然完全下不了手&#xff0c;本文就以IEEE的模板&#xff0c;从头往下讲咯~ 第…

【代码随想录】【算法训练营】【第44天】 [322]零钱兑换 [279]完全平方数 [139]单词拆分

前言 思路及算法思维&#xff0c;指路 代码随想录。 题目来自 LeetCode。 day 44&#xff0c;周四&#xff0c;坚持不住了~ 题目详情 [322] 零钱兑换 题目描述 322 零钱兑换 解题思路 前提&#xff1a; 思路&#xff1a; 重点&#xff1a; 代码实现 C语言 [279] 完全…

7亿中国男人,今年夏天都在穿什么?

文丨郭梦仪 北京气温已经逼近38度&#xff0c;注重防晒的人群中这次多了男人的身影。 程序员宇宙中心&#xff0c;清河万象汇西区&#xff0c;小米su7吸引众多男士前来观摩&#xff0c;和对面蕉下门店里的“防晒衣大军”恰好呼应上了。 北京清河万象汇的防晒衣专卖店 夏日将…

Studying-代码随想录训练营day15| 222.完全二叉树的节点个数、110.平衡二叉树、257.二叉树的所有路径、404.左叶子之和

第十五天&#xff0c;二叉树part03&#x1f4aa;&#xff0c;编程语言&#xff1a;C 目录 257.完全二叉树的节点个数 110.平衡二叉树 257.二叉树的所有路径 404.左叶子之和 总结 257.完全二叉树的节点个数 文档讲解&#xff1a;代码随想录完全二叉树的节点个数 视频讲解…

118 杨辉三角

题目 给定一个非负整数 numRows&#xff0c;生成「杨辉三角」的前 numRows 行。 在「杨辉三角」中&#xff0c;每个数是它左上方和右上方的数的和。 示例 输入: numRows 5 输出: [[1],[1,1],[1,2,1],[1,3,3,1],[1,4,6,4,1]] 解析 就是模拟法&#xff0c;没有什么特殊的…

从仓位出发 略谈在线伦敦银交易的技巧

现在我们做伦敦银投资&#xff0c;都是通过线上完成的。在线做伦敦金交易的意思&#xff0c;就是投资者通过网络&#xff0c;在PC端或者移动端去利用交易软件完成交易&#xff0c;那么在线伦敦银交易有什么技巧呢&#xff1f;下面我们就从仓位的角度来讨论。 投资者入场后所持有…

Vulhub——Log4j、solr

文章目录 一、Log4j1.1 Apache Log4j2 lookup JNDI 注入漏洞&#xff08;CVE-2021-44228&#xff09;1.2 Apache Log4j Server 反序列化命令执行漏洞&#xff08;CVE-2017-5645&#xff09; 二、Solr2.1 Apache Solr 远程命令执行漏洞&#xff08;CVE-2017-12629&#xff09;2.…

Java 笔记:常见正则使用

文章目录 Java 笔记&#xff1a;常见正则使用正则简介常用匹配年月日的时间匹配手机号码校验 参考文章 Java 笔记&#xff1a;常见正则使用 正则简介 正则表达式定义了字符串的模式。 正则表达式可以用来搜索、编辑或处理文本。 正则表达式并不仅限于某一种语言&#xff0c;但…

Freemaker 模板

背景 发送邮件&#xff0c;正文利用freemaker完成 <dependency><groupId>org.springframework.boot</groupId><artifactId>spring-boot-starter-freemarker</artifactId> </dependency>Autowired private Configuration configurer;GetMap…

Set集合系列——Set、HashSet、LinkedHashset、TreeSet

Set系列的公共特点&#xff1a;无重复、无索引&#xff0c;不可用普通for循环&#xff0c;API和Collection重复 HashSet&#xff1a;采取哈希表存取数据 哈希表组成&#xff1f; JDk8之前&#xff1a;数组链表&#xff0c; JDK8以后&#xff1a;数组链表红黑树 哈希值&#…

【Python】使用matplotlib绘制图形(曲线图、条形图、饼图等)

文章目录 一、什么是matplotlib二、matplotlib 支持的图形三、如何使用matplotlib1. 安装matplotlib2. 导入matplotlib.pyplot3. 准备数据4. 绘制图形5. 定制图形6. 显示或保存图形7. &#xff08;可选&#xff09;使用subplots创建多个子图注意事项&#xff1a; 四、常见图形使…

CCF推荐会议必投攻略:这些顶级会议投完直通录取大门

CCF推荐会议必投攻略&#xff1a;这些顶级会议投完直通录取大门&#xff01; 会议之眼 快讯 CCF介绍 CCF&#xff08;China Computer Federation&#xff09;即中国计算机学会&#xff0c;前身是中国电子学会计算机专业委员会&#xff0c;成立于1962年。这是由从事计算机及相…

idea2023开发插件入门

idea2023开发插件入门 创建工程 通过 idea plugin 来创建工程 修改 开发语言 默认创建的工程是用scala开发的&#xff0c;但是我不会&#xff0c;就会java,所以改成java创建 build.gradle.kt 为 build.gradlesettings.gradle.kt 为 settings.gradle build.gradle修改为以…

食品安全无小事:EasyCVR+AI技术助力食品加工厂管理透明化,构建食品安全防线

一、背景需求 近期有新闻记者曝光某省禽类屠宰加工厂脏乱差问题严重&#xff0c;工人脚踩鹅肠鸭肠混杂洗地水、烟头随手扔进鸭肠筐、污水捞出死鸭再上生产线…卫生情况十分堪忧。食品卫生安全频频出现负面新闻&#xff0c;如何实现源头治理&#xff1f;如何将各类食品安全风险隐…

聊聊 oracle varchar2 字段的gbk/utf8编码格式和字段长度问题

聊聊 oracle varchar2 字段的gbk/utf8编码格式和字段长度问题 1 问题现象 最近在排查某客户现场的数据同步作业报错问题时&#xff0c;发现了部分 ORACLE 表的 varchar2 字段&#xff0c;因为上游 ORACLE数据库采用 GBK 编码格式&#xff0c;而下游 ORACLE 数据库采用UTF8 编…

开发一个软件自动运行工具不可缺少的源代码分享!

在软件开发领域&#xff0c;自动运行工具扮演着至关重要的角色&#xff0c;它们能够简化软件部署、提升运行效率&#xff0c;并在很大程度上降低人为操作失误的可能性。 而一个高效的自动运行工具的背后&#xff0c;往往是经过精心设计与实现的源代码在默默支撑&#xff0c;本…

html做一个画柱形图的软件

你可以使用 HTML、CSS 和 JavaScript 创建一个简单的柱形图绘制软件。为了方便起见&#xff0c;我们可以使用一个流行的 JavaScript 图表库&#xff0c;比如 Chart.js&#xff0c;它能够简化创建和操作图表的过程。 以下是一个完整的示例&#xff0c;展示如何使用 HTML 和 Cha…

代码随想录-Day36

452. 用最少数量的箭引爆气球 有一些球形气球贴在一堵用 XY 平面表示的墙面上。墙面上的气球记录在整数数组 points &#xff0c;其中points[i] [xstart, xend] 表示水平直径在 xstart 和 xend之间的气球。你不知道气球的确切 y 坐标。 一支弓箭可以沿着 x 轴从不同点 完全垂…