目录
介绍
主要步骤
代码
__init__.py
min_bounding_rect.py
min_rect.py
qhull_2d.py
结果
介绍
最小有向包围盒算法广泛应用于多个领域，包括：
- 计算几何：用于分析点集的边界特征。
- 图形学：用于碰撞检测和物体包围。
- 数据科学：在聚类分析、异常值检测等场景中用于计算数据的紧致边界。
本文提出了一种基于凸包和最小矩形算法的二维最小有向包围盒计算方法，并展示了如何通过代码实现这一算法。
主要步骤
- 凸包计算：使用
qhull2D方法计算输入点集的凸包。这一步可以确保我们找到点集的外部边界，以此作为最小矩形的候选边界。 - 最小矩形计算：利用
minBoundingRect算法计算最小有向包围盒。算法返回矩形的宽度、高度、角点坐标和面积。
程序结构如下：
调用入口：
输入为numpy格式的n*2数组。
代码
__init__.py
from .min_rect import *
from .qhull_2d import *
from .min_bounding_rect import *
min_bounding_rect.py
import numpy as np
import sysdef minBoundingRect(hull_points_2d):# Compute edges (x2-x1, y2-y1)edges = np.zeros((len(hull_points_2d)-1, 2)) # empty 2 column arrayfor i in range(len(edges)):edge_x = hull_points_2d[i+1, 0] - hull_points_2d[i, 0]edge_y = hull_points_2d[i+1, 1] - hull_points_2d[i, 1]edges[i] = [edge_x, edge_y]# Calculate edge angles atan2(y/x)edge_angles = np.zeros(len(edges)) # empty 1 column arrayfor i in range(len(edge_angles)):edge_angles[i] = np.arctan2(edges[i, 1], edges[i, 0])# Check for angles in 1st quadrantedge_angles = np.abs(edge_angles % (np.pi/2)) # want strictly positive answers# Remove duplicate anglesedge_angles = np.unique(edge_angles)# Test each angle to find bounding box with smallest areamin_bbox = (0, sys.maxsize, 0, 0, 0, 0, 0, 0) # rot_angle, area, width, height, min_x, max_x, min_y, max_y# print("Testing", len(edge_angles), "possible rotations for bounding box... \n")for i in range(len(edge_angles)):# Create rotation matrix to shift points to baselineR = np.array([[np.cos(edge_angles[i]), np.cos(edge_angles[i]-(np.pi/2))],[np.cos(edge_angles[i]+(np.pi/2)), np.cos(edge_angles[i])]])# Apply this rotation to convex hull pointsrot_points = np.dot(R, np.transpose(hull_points_2d)) # 2x2 * 2xn# Find min/max x,y pointsmin_x = np.nanmin(rot_points[0], axis=0)max_x = np.nanmax(rot_points[0], axis=0)min_y = np.nanmin(rot_points[1], axis=0)max_y = np.nanmax(rot_points[1], axis=0)# Calculate height/width/area of this bounding rectanglewidth = max_x - min_xheight = max_y - min_yarea = width * height# Store the smallest rect found firstif area < min_bbox[1]:min_bbox = (edge_angles[i], area, width, height, min_x, max_x, min_y, max_y)# Re-create rotation matrix for smallest rectangle = min_bbox[0]R = np.array([[np.cos(angle), np.cos(angle-(np.pi/2))],[np.cos(angle+(np.pi/2)), np.cos(angle)]])# Project convex hull points onto rotated frameproj_points = np.dot(R, np.transpose(hull_points_2d)) # 2x2 * 2xn# min/max x,y points are against baselinemin_x = min_bbox[4]max_x = min_bbox[5]min_y = min_bbox[6]max_y = min_bbox[7]# Calculate center point and project onto rotated frame# center_x = (min_x + max_x) / 2# center_y = (min_y + max_y) / 2# center_point = np.dot([center_x, center_y], R)# Calculate corner points and project onto rotated framecorner_points = np.zeros((4, 2)) # empty 2 column arraycorner_points[0] = np.dot([max_x, min_y], R)corner_points[1] = np.dot([min_x, min_y], R)corner_points[2] = np.dot([min_x, max_y], R)corner_points[3] = np.dot([max_x, max_y], R)return min_bbox[2], min_bbox[3], corner_points,areamin_rect.py
import numpy as np
import matplotlib.pyplot as plt
from .qhull_2d import qhull2D
from .min_bounding_rect import minBoundingRectdef min_area_rect(xy_points):# Find convex hullhull_points = qhull2D(xy_points)# Reverse order of points, to match output from other qhull implementationshull_points = hull_points[::-1]# print('Convex hull points: \n', hull_points, "\n")# Find minimum area bounding rectanglewidth, height, corner_points,area = minBoundingRect(hull_points)# # Visualizationplt.figure()plt.plot(xy_points[:, 0], xy_points[:, 1], 'o', label='Points')plt.plot(hull_points[:, 0], hull_points[:, 1], 'r--', lw=2, label='Convex Hull')# Plot the bounding boxbox = np.vstack([corner_points, corner_points[0]]) # Close the box by repeating the first pointplt.plot(box[:, 0], box[:, 1], 'g-', lw=2, label='Min Bounding Box')plt.legend()plt.xlabel('X')plt.ylabel('Y')plt.title('Convex Hull and Minimum Area Bounding Box')plt.grid(True)plt.axis('equal')plt.show()return width, height,corner_points,areaqhull_2d.py
from __future__ import division
from numpy import *def link(a, b):return concatenate((a, b[1:]))def edge(a, b):return concatenate(([a], [b]))def qhull2D(sample):def dome(sample, base, depth=0, max_depth=1000):h, t = basedists = dot(sample - h, dot(((0, -1), (1, 0)), (t - h)))if len(dists) == 0:return base# Handle cases where all distances are very smallif all(abs(dists) < 1e-10):return baseouter = sample[dists > 0]# print("dists:",dists,"outer:",outer)# Handle empty outer caseif len(outer) == 0:return basepivot_index = argmax(dists)pivot = sample[pivot_index]# Ensure depth does not exceed maximum allowedif depth > max_depth:return base# Recursive caseleft_dome = dome(outer, edge(h, pivot), depth + 1, max_depth)right_dome = dome(outer, edge(pivot, t), depth + 1, max_depth)return link(left_dome, right_dome)if len(sample) > 2:axis = sample[:, 0]base = take(sample, [argmin(axis), argmax(axis)], axis=0)left_dome = dome(sample, base)right_dome = dome(sample, base[::-1])return link(left_dome, right_dome)else:return sample
结果
