|
12 | 12 | __all__ = [ |
13 | 13 | 'bounding_box', |
14 | 14 | 'bounding_box_xy', |
15 | | - 'oriented_bounding_box_numpy', |
16 | | - 'oriented_bounding_box_xy_numpy', |
17 | 15 | ] |
18 | 16 |
|
19 | 17 |
|
@@ -89,245 +87,6 @@ def bounding_box_xy(points): |
89 | 87 | (min_x, max_y, 0.0)] |
90 | 88 |
|
91 | 89 |
|
92 | | -def oriented_bounding_box_numpy(points): |
93 | | - """Compute the oriented minimum bounding box of a set of points in 3D space. |
94 | | -
|
95 | | - Notes |
96 | | - ----- |
97 | | - The implementation is based on the convex hull of the points. |
98 | | -
|
99 | | - Parameters |
100 | | - ---------- |
101 | | - points : list |
102 | | - XYZ coordinates of the points. |
103 | | -
|
104 | | - Returns |
105 | | - ------- |
106 | | - list |
107 | | - The XYZ coordinates of the corners of the bounding box. |
108 | | -
|
109 | | - Examples |
110 | | - -------- |
111 | | - .. plot:: |
112 | | - :include-source: |
113 | | -
|
114 | | - from numpy.random import randint |
115 | | - from numpy.random import rand |
116 | | -
|
117 | | - import matplotlib.pyplot as plt |
118 | | -
|
119 | | - from compas.plotters import Bounds |
120 | | - from compas.plotters import Cloud3D |
121 | | - from compas.plotters import Box |
122 | | - from compas.plotters import create_axes_3d |
123 | | - from compas.geometry import rotation_matrix |
124 | | - from compas.geometry import transform |
125 | | - from compas.geometry import oriented_bounding_box_numpy |
126 | | -
|
127 | | - clouds = [] |
128 | | -
|
129 | | - for i in range(8): |
130 | | - a = randint(1, high=8) * 10 * 3.14159 / 180 |
131 | | - d = [1, 1, 1] |
132 | | -
|
133 | | - cloud = rand(100, 3) |
134 | | -
|
135 | | - if i in (1, 2, 5, 6): |
136 | | - cloud[:, 0] *= - 10.0 |
137 | | - cloud[:, 0] -= 3.0 |
138 | | - d[0] = -1 |
139 | | - else: |
140 | | - cloud[:, 0] *= 10.0 |
141 | | - cloud[:, 0] += 3.0 |
142 | | -
|
143 | | - if i in (2, 3, 6, 7): |
144 | | - cloud[:, 1] *= - 3.0 |
145 | | - cloud[:, 1] -= 3.0 |
146 | | - d[1] = -1 |
147 | | - else: |
148 | | - cloud[:, 1] *= 3.0 |
149 | | - cloud[:, 1] += 3.0 |
150 | | -
|
151 | | - if i in (4, 5, 6, 7): |
152 | | - cloud[:, 2] *= - 6.0 |
153 | | - cloud[:, 2] -= 3.0 |
154 | | - d[2] = -1 |
155 | | - else: |
156 | | - cloud[:, 2] *= 6.0 |
157 | | - cloud[:, 2] += 3.0 |
158 | | -
|
159 | | - R = rotation_matrix(a, d) |
160 | | - cloud[:] = transform(cloud, R) |
161 | | -
|
162 | | - clouds.append(cloud.tolist()) |
163 | | -
|
164 | | - axes = create_axes_3d() |
165 | | -
|
166 | | - bounds = Bounds([point for points in clouds for point in points]) |
167 | | - bounds.plot(axes) |
168 | | -
|
169 | | - for cloud in clouds: |
170 | | - bbox = oriented_bounding_box_numpy(cloud) |
171 | | -
|
172 | | - Cloud3D(cloud).plot(axes) |
173 | | - Box(bbox[1]).plot(axes) |
174 | | -
|
175 | | - plt.show() |
176 | | -
|
177 | | - """ |
178 | | - from numpy import asarray |
179 | | - from numpy import argmax |
180 | | - from numpy import argmin |
181 | | - from numpy import ptp |
182 | | - from scipy.spatial import ConvexHull |
183 | | - |
184 | | - from compas.geometry import local_axes |
185 | | - from compas.geometry import local_coords_numpy |
186 | | - from compas.geometry import global_coords_numpy |
187 | | - |
188 | | - points = asarray(points) |
189 | | - n, dim = points.shape |
190 | | - |
191 | | - assert 2 < dim, "The point coordinates should be at least 3D: %i" % dim |
192 | | - |
193 | | - points = points[:, :3] |
194 | | - |
195 | | - hull = ConvexHull(points) |
196 | | - volume = None |
197 | | - bbox = [] |
198 | | - |
199 | | - # this can be vectorised! |
200 | | - for simplex in hull.simplices: |
201 | | - abc = points[simplex] |
202 | | - uvw = local_axes(abc[0], abc[1], abc[2]) |
203 | | - xyz = points[hull.vertices] |
204 | | - rst = local_coords_numpy(abc[0], uvw, xyz) |
205 | | - dr, ds, dt = ptp(rst, axis=0) |
206 | | - v = dr * ds * dt |
207 | | - |
208 | | - if volume is None or v < volume: |
209 | | - rmin, smin, tmin = argmin(rst, axis=0) |
210 | | - rmax, smax, tmax = argmax(rst, axis=0) |
211 | | - bbox = [ |
212 | | - [rst[rmin, 0], rst[smin, 1], rst[tmin, 2]], |
213 | | - [rst[rmax, 0], rst[smin, 1], rst[tmin, 2]], |
214 | | - [rst[rmax, 0], rst[smax, 1], rst[tmin, 2]], |
215 | | - [rst[rmin, 0], rst[smax, 1], rst[tmin, 2]], |
216 | | - [rst[rmin, 0], rst[smin, 1], rst[tmax, 2]], |
217 | | - [rst[rmax, 0], rst[smin, 1], rst[tmax, 2]], |
218 | | - [rst[rmax, 0], rst[smax, 1], rst[tmax, 2]], |
219 | | - [rst[rmin, 0], rst[smax, 1], rst[tmax, 2]], |
220 | | - ] |
221 | | - bbox = global_coords_numpy(abc[0], uvw, bbox) |
222 | | - volume = v |
223 | | - |
224 | | - return hull, bbox, volume |
225 | | - |
226 | | - |
227 | | -def oriented_bounding_box_xy_numpy(points): |
228 | | - """Compute the oriented minimum bounding box of set of points in the XY plane. |
229 | | -
|
230 | | - Notes |
231 | | - ----- |
232 | | - The *oriented (minimum) bounding box* (OBB) is computed using the following |
233 | | - procedure: |
234 | | -
|
235 | | - 1. Compute the convex hull of the points. |
236 | | - 2. For each of the edges on the hull: |
237 | | - 1. Compute the s-axis as the unit vector in the direction of the edge |
238 | | - 2. Compute the othorgonal t-axis. |
239 | | - 3. Use the start point of the edge as origin. |
240 | | - 4. Compute the spread of the points along the s-axis. |
241 | | - (dot product of the point vecor in local coordinates and the s-axis) |
242 | | - 5. Compute the spread along the t-axis. |
243 | | - 6. Determine the side of s on which the points are. |
244 | | - 7. Compute and store the corners of the bbox and its area. |
245 | | - 3. Select the box with the smallest area. |
246 | | -
|
247 | | - Parameters |
248 | | - ---------- |
249 | | - points : list |
250 | | - XY(Z) coordinates of the points. |
251 | | -
|
252 | | - Returns |
253 | | - ------- |
254 | | - list |
255 | | - The coordinates of the corners of the bounding box. |
256 | | -
|
257 | | - Examples |
258 | | - -------- |
259 | | - .. code-block:: python |
260 | | -
|
261 | | - # |
262 | | -
|
263 | | - """ |
264 | | - from numpy import array |
265 | | - from numpy import asarray |
266 | | - from numpy import argmax |
267 | | - from numpy import argmin |
268 | | - from numpy import dot |
269 | | - |
270 | | - from scipy.spatial import ConvexHull |
271 | | - |
272 | | - points = asarray(points) |
273 | | - n, dim = points.shape |
274 | | - |
275 | | - assert 1 < dim, "The point coordinates should be at least 2D: %i" % dim |
276 | | - |
277 | | - points = points[:, :2] |
278 | | - |
279 | | - hull = ConvexHull(points) |
280 | | - xy_hull = points[hull.vertices].reshape((-1, 2)) |
281 | | - |
282 | | - boxes = [] |
283 | | - m = sum(xy_hull, axis=0) / n |
284 | | - |
285 | | - for simplex in hull.simplices: |
286 | | - p0 = points[simplex[0]] |
287 | | - p1 = points[simplex[1]] |
288 | | - |
289 | | - # s direction |
290 | | - s = p1 - p0 |
291 | | - sl = sum(s ** 2) ** 0.5 |
292 | | - su = s / sl |
293 | | - vn = xy_hull - p0 |
294 | | - sc = (sum(vn * s, axis=1) / sl).reshape((-1, 1)) |
295 | | - scmax = argmax(sc) |
296 | | - scmin = argmin(sc) |
297 | | - |
298 | | - # box corners |
299 | | - b0 = p0 + sc[scmin] * su |
300 | | - b1 = p0 + sc[scmax] * su |
301 | | - |
302 | | - # t direction |
303 | | - t = array([-s[1], s[0]]) |
304 | | - tl = sum(t ** 2) ** 0.5 |
305 | | - tu = t / tl |
306 | | - vn = xy_hull - p0 |
307 | | - tc = (sum(vn * t, axis=1) / tl).reshape((-1, 1)) |
308 | | - tcmax = argmax(tc) |
309 | | - tcmin = argmin(tc) |
310 | | - |
311 | | - # area |
312 | | - w = sc[scmax] - sc[scmin] |
313 | | - h = tc[tcmax] - tc[tcmin] |
314 | | - a = w * h |
315 | | - |
316 | | - # box corners |
317 | | - if dot(t, m - p0) < 0: |
318 | | - b3 = b0 - h * tu |
319 | | - b2 = b1 - h * tu |
320 | | - else: |
321 | | - b3 = b0 + h * tu |
322 | | - b2 = b1 + h * tu |
323 | | - |
324 | | - # box |
325 | | - boxes.append([[b0, b1, b2, b3], a[0]]) |
326 | | - |
327 | | - # return the box with the smallest area |
328 | | - return min(boxes, key=lambda b: b[1]) |
329 | | - |
330 | | - |
331 | 90 | # ============================================================================== |
332 | 91 | # Main |
333 | 92 | # ============================================================================== |
|
0 commit comments