计算多个凸二维多边形交集的建议
Suggestions to Compute the Intersetions of Multiple Convex 2D Polygons
我写这个问题是为了寻找可以快速计算 N
二维多边形(投影凸多面体的凸包)和 M
二维多边形,通常为 N >> M
。 N
可能在数量级或至少 1M 多边形和 N
在数量级 50k。我已经搜索了一段时间,但我不断得出如下所示的相同答案。
使用 boost 和循环
- 计算多面体的投影(不是瓶颈)
- 计算所述多面体的凸包(瓶颈)
- 计算投影多面体与现有二维多边形的交集(主要瓶颈)。
此循环重复 NK
次,其中通常 K << M
,K
是与单个投影多面体相交的二维多边形的平均数量。这样做是为了减少计算量。
这个问题是,如果我有 N=262144
和 M=19456
,它需要大约 129 秒(当多面体多线程时),并且必须完成大约 300 次。理想情况下,我想将上述尺寸的计算时间减少到大约 1 秒,所以我想知道是否有人可以帮助指出一些可以提高效率的软件或文献。
[编辑]
@sehe 的请求我发布了代码中最相关的部分。我还没有编译它,所以这只是为了了解要点……这段代码假设有体素和像素,但形状可以是任何东西。网格中点的顺序可以任意,但点在网格中的索引是相同的。
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/point.hpp>
#include <boost/geometry/geometries/ring.hpp>
const std::size_t Dimension = 2;
typedef boost::geometry::model::point<float, Dimension, boost::geometry::cs::cartesian> point_2d;
typedef boost::geometry::model::polygon<point_2d, false /* is cw */, true /* closed */> polygon_2d;
typedef boost::geometry::model::box<point_2d> box_2d;
std::vector<float> getOverlaps(std::vector<float> & projected_grid_vx, // projected voxels
std::vector<float> & pixel_grid_vx, // pixels
std::vector<int> & projected_grid_m, // number of voxels in each dimension
std::vector<int> & pixel_grid_m, // number of pixels in each dimension
std::vector<float> & pixel_grid_omega, // size of the pixel grid in cm
int projected_grid_size, // total number of voxels
int pixel_grid_size) { // total number of pixels
std::vector<float> overlaps(projected_grid_size * pixel_grid_size);
std::vector<float> h(pixel_grid_m.size());
for(int d=0; d < pixel_grid_m.size(); d++) {
h[d] = (pixel_grid_omega[2*d+1] - pixel_grid_omega[2*d]) / pixel_grid_m[d];
}
for(int i=0; i < projected_grid_size; i++){
std::vector<float> point_indices(8);
point_indices[0] = i;
point_indices[1] = i + 1;
point_indices[2] = i + projected_grid_m[0];
point_indices[3] = i + projected_grid_m[0] + 1;
point_indices[4] = i + projected_grid_m[0] * projected_grid_m[1];
point_indices[5] = i + projected_grid_m[0] * projected_grid_m[1] + 1;
point_indices[6] = i + (projected_grid_m[1] + 1) * projected_grid_m[0];
point_indices[7] = i + (projected_grid_m[1] + 1) * projected_grid_m[0] + 1;
std::vector<float> vx_corners(8 * projected_grid_m.size());
for(int vn = 0; vn < 8; vn++) {
for(int d = 0; d < projected_grid_m.size(); d++) {
vx_corners[vn + d * 8] = projected_grid_vx[point_indices[vn] + d * projeted_grid_size];
}
}
polygon_2d proj_voxel;
for(int vn = 0; vn < 8; vn++) {
point_2d poly_pt(vx_corners[2 * vn], vx_corners[2 * vn + 1]);
boost::geometry::append(proj_voxel, poly_pt);
}
boost::geometry::correct(proj_voxel);
polygon_2d proj_voxel_hull;
boost::geometry::convex_hull(proj_voxel, proj_voxel_hull);
box_2d bb_proj_vox;
boost::geometry::envelope(proj_voxel_hull, bb_proj_vox);
point_2d min_pt = bb_proj_vox.min_corner();
point_2d max_pt = bb_proj_vox.max_corner();
// then get min and max indices of intersecting bins
std::vector<float> min_idx(projected_grid_m.size() - 1),
max_idx(projected_grid_m.size() - 1);
// compute min and max indices of incidence on the pixel grid
// this is easy assuming you have a regular grid of pixels
min_idx[0] = std::min( (float) std::max( std::floor((min_pt.get<0>() - pixel_grid_omega[0]) / h[0] - 0.5 ), 0.), pixel_grid_m[0]-1);
min_idx[1] = std::min( (float) std::max( std::floor((min_pt.get<1>() - pixel_grid_omega[2]) / h[1] - 0.5 ), 0.), pixel_grid_m[1]-1);
max_idx[0] = std::min( (float) std::max( std::floor((max_pt.get<0>() - pixel_grid_omega[0]) / h[0] + 0.5 ), 0.), pixel_grid__m[0]-1);
max_idx[1] = std::min( (float) std::max( std::floor((max_pt.get<1>() - pixel_grid_omega[2]) / h[1] + 0.5 ), 0.), pixel_grid_m[1]-1);
// iterate only over pixels which intersect the projected voxel
for(int iy = min_idx[1]; iy <= max_idx[1]; iy++) {
for(int ix = min_idx[0]; ix <= max_idx[0]; ix++) {
int idx = ix + iy * pixel_grid_size[0]; // `first' index of pixel corner point
polygon_2d pix_poly;
for(int pn = 0; pn < 4; pn++) {
point_2d pix_corner_pt(
pixel_grid_vx[idx + pn % 2 + (pn / 2) * pixel_grid_m[0]],
pixel_grid_vx[idx + pn % 2 + (pn / 2) * pixel_grid_m[0] + pixel_grid_size]
);
boost::geometry::append(pix_poly, pix_corner_pt);
}
boost::geometry::correct( pix_poly );
//make this into a convex hull since the order of the point may be any
polygon_2d pix_hull;
boost::geometry::convex_hull(pix_poly, pix_hull);
// on to perform intersection
std::vector<polygon_2d> vox_pix_ints;
polygon_2d vox_pix_int;
try {
boost::geometry::intersection(proj_voxel_hull, pix_hull, vox_pix_ints);
} catch ( std::exception e ) {
// skip since these may coincide at a point or line
continue;
}
// both are convex so only one intersection expected
vox_pix_int = vox_pix_ints[0];
overlaps[i + idx * projected_grid_size] = boost::geometry::area(vox_pix_int);
}
} // end intersection for
} //end projected_voxel for
return overlaps;
}
嗯,问题似乎与 "collision-detection"i 游戏引擎类似。或者 "potentially visible sets"。
虽然我对当前的最新技术知之甚少,但我记得优化是将对象封闭在球体中,因为检查球体(或 2D 中的圆)之间的重叠非常便宜。
为了加快碰撞检查,通常将对象放入搜索结构中(例如球体树(二维情况下的圆树))。基本上将 space 组织成层次结构,以便快速查询重叠。
基本上我的建议可以归结为:尝试查看游戏引擎中的碰撞检测算法。
假设
我假设你的意思是 "intersections" 而不是交集。此外,M
和 N
中的大多数单个多边形同时重叠并不是预期的用例。如果这个假设成立,那么:
回答
使用 2D 游戏引擎完成此操作的方法是使用场景图,其中每个对象都有一个边界框。然后根据边界框确定的位置将所有多边形放入四叉树中的一个节点。那么任务就变成了并行,因为每个节点都可以单独处理交集。
这是四叉树的 wiki:
在 3D 中可以使用八叉树。
它实际上什至不必是八叉树。您可以使用任何 space 分区获得相同的结果。你可以找到多边形的最大分离(我们称之为 S
)。并创建说 S/10
space 分区。然后您将有 10 个单独的 space 并行执行。它不仅是并发的,而且不再是 M * N
时间,因为不是每个多边形都必须与其他所有多边形进行比较。
您可以创建多边形与边界框的比率:
这可以在计算上完成一次,以获得平均多边形区域与 BB 比率 R
常数。
或者您可以使用以 BB 为界的圆来处理几何因为您只使用投影多面体:
R = 0.0;
count = 0;
for (each poly) {
count++;
R += polyArea / itsBoundingBoxArea;
}
R = R/count;
然后计算边界框交集的总和。
Sbb = 0.0;
for (box1, box2 where box1.isIntersecting(box2)) {
Sbb += box1.intersect(box2);
}
然后:
Approximation = R * Sbb
如果允许凹多边形,所有这些都将不起作用。因为凹多边形可以占据不到其边界框的 1%。您仍然需要找到凸包。
或者,如果您可以比其外壳更快地找到多边形面积,则可以使用实际计算的平均多边形面积。这也会给你一个不错的近似值,同时避免多边形相交和环绕。
我写这个问题是为了寻找可以快速计算 N
二维多边形(投影凸多面体的凸包)和 M
二维多边形,通常为 N >> M
。 N
可能在数量级或至少 1M 多边形和 N
在数量级 50k。我已经搜索了一段时间,但我不断得出如下所示的相同答案。
使用 boost 和循环
- 计算多面体的投影(不是瓶颈)
- 计算所述多面体的凸包(瓶颈)
- 计算投影多面体与现有二维多边形的交集(主要瓶颈)。
此循环重复 NK
次,其中通常 K << M
,K
是与单个投影多面体相交的二维多边形的平均数量。这样做是为了减少计算量。
这个问题是,如果我有 N=262144
和 M=19456
,它需要大约 129 秒(当多面体多线程时),并且必须完成大约 300 次。理想情况下,我想将上述尺寸的计算时间减少到大约 1 秒,所以我想知道是否有人可以帮助指出一些可以提高效率的软件或文献。
[编辑]
@sehe 的请求我发布了代码中最相关的部分。我还没有编译它,所以这只是为了了解要点……这段代码假设有体素和像素,但形状可以是任何东西。网格中点的顺序可以任意,但点在网格中的索引是相同的。
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/point.hpp>
#include <boost/geometry/geometries/ring.hpp>
const std::size_t Dimension = 2;
typedef boost::geometry::model::point<float, Dimension, boost::geometry::cs::cartesian> point_2d;
typedef boost::geometry::model::polygon<point_2d, false /* is cw */, true /* closed */> polygon_2d;
typedef boost::geometry::model::box<point_2d> box_2d;
std::vector<float> getOverlaps(std::vector<float> & projected_grid_vx, // projected voxels
std::vector<float> & pixel_grid_vx, // pixels
std::vector<int> & projected_grid_m, // number of voxels in each dimension
std::vector<int> & pixel_grid_m, // number of pixels in each dimension
std::vector<float> & pixel_grid_omega, // size of the pixel grid in cm
int projected_grid_size, // total number of voxels
int pixel_grid_size) { // total number of pixels
std::vector<float> overlaps(projected_grid_size * pixel_grid_size);
std::vector<float> h(pixel_grid_m.size());
for(int d=0; d < pixel_grid_m.size(); d++) {
h[d] = (pixel_grid_omega[2*d+1] - pixel_grid_omega[2*d]) / pixel_grid_m[d];
}
for(int i=0; i < projected_grid_size; i++){
std::vector<float> point_indices(8);
point_indices[0] = i;
point_indices[1] = i + 1;
point_indices[2] = i + projected_grid_m[0];
point_indices[3] = i + projected_grid_m[0] + 1;
point_indices[4] = i + projected_grid_m[0] * projected_grid_m[1];
point_indices[5] = i + projected_grid_m[0] * projected_grid_m[1] + 1;
point_indices[6] = i + (projected_grid_m[1] + 1) * projected_grid_m[0];
point_indices[7] = i + (projected_grid_m[1] + 1) * projected_grid_m[0] + 1;
std::vector<float> vx_corners(8 * projected_grid_m.size());
for(int vn = 0; vn < 8; vn++) {
for(int d = 0; d < projected_grid_m.size(); d++) {
vx_corners[vn + d * 8] = projected_grid_vx[point_indices[vn] + d * projeted_grid_size];
}
}
polygon_2d proj_voxel;
for(int vn = 0; vn < 8; vn++) {
point_2d poly_pt(vx_corners[2 * vn], vx_corners[2 * vn + 1]);
boost::geometry::append(proj_voxel, poly_pt);
}
boost::geometry::correct(proj_voxel);
polygon_2d proj_voxel_hull;
boost::geometry::convex_hull(proj_voxel, proj_voxel_hull);
box_2d bb_proj_vox;
boost::geometry::envelope(proj_voxel_hull, bb_proj_vox);
point_2d min_pt = bb_proj_vox.min_corner();
point_2d max_pt = bb_proj_vox.max_corner();
// then get min and max indices of intersecting bins
std::vector<float> min_idx(projected_grid_m.size() - 1),
max_idx(projected_grid_m.size() - 1);
// compute min and max indices of incidence on the pixel grid
// this is easy assuming you have a regular grid of pixels
min_idx[0] = std::min( (float) std::max( std::floor((min_pt.get<0>() - pixel_grid_omega[0]) / h[0] - 0.5 ), 0.), pixel_grid_m[0]-1);
min_idx[1] = std::min( (float) std::max( std::floor((min_pt.get<1>() - pixel_grid_omega[2]) / h[1] - 0.5 ), 0.), pixel_grid_m[1]-1);
max_idx[0] = std::min( (float) std::max( std::floor((max_pt.get<0>() - pixel_grid_omega[0]) / h[0] + 0.5 ), 0.), pixel_grid__m[0]-1);
max_idx[1] = std::min( (float) std::max( std::floor((max_pt.get<1>() - pixel_grid_omega[2]) / h[1] + 0.5 ), 0.), pixel_grid_m[1]-1);
// iterate only over pixels which intersect the projected voxel
for(int iy = min_idx[1]; iy <= max_idx[1]; iy++) {
for(int ix = min_idx[0]; ix <= max_idx[0]; ix++) {
int idx = ix + iy * pixel_grid_size[0]; // `first' index of pixel corner point
polygon_2d pix_poly;
for(int pn = 0; pn < 4; pn++) {
point_2d pix_corner_pt(
pixel_grid_vx[idx + pn % 2 + (pn / 2) * pixel_grid_m[0]],
pixel_grid_vx[idx + pn % 2 + (pn / 2) * pixel_grid_m[0] + pixel_grid_size]
);
boost::geometry::append(pix_poly, pix_corner_pt);
}
boost::geometry::correct( pix_poly );
//make this into a convex hull since the order of the point may be any
polygon_2d pix_hull;
boost::geometry::convex_hull(pix_poly, pix_hull);
// on to perform intersection
std::vector<polygon_2d> vox_pix_ints;
polygon_2d vox_pix_int;
try {
boost::geometry::intersection(proj_voxel_hull, pix_hull, vox_pix_ints);
} catch ( std::exception e ) {
// skip since these may coincide at a point or line
continue;
}
// both are convex so only one intersection expected
vox_pix_int = vox_pix_ints[0];
overlaps[i + idx * projected_grid_size] = boost::geometry::area(vox_pix_int);
}
} // end intersection for
} //end projected_voxel for
return overlaps;
}
嗯,问题似乎与 "collision-detection"i 游戏引擎类似。或者 "potentially visible sets"。
虽然我对当前的最新技术知之甚少,但我记得优化是将对象封闭在球体中,因为检查球体(或 2D 中的圆)之间的重叠非常便宜。 为了加快碰撞检查,通常将对象放入搜索结构中(例如球体树(二维情况下的圆树))。基本上将 space 组织成层次结构,以便快速查询重叠。
基本上我的建议可以归结为:尝试查看游戏引擎中的碰撞检测算法。
假设
我假设你的意思是 "intersections" 而不是交集。此外,M
和 N
中的大多数单个多边形同时重叠并不是预期的用例。如果这个假设成立,那么:
回答
使用 2D 游戏引擎完成此操作的方法是使用场景图,其中每个对象都有一个边界框。然后根据边界框确定的位置将所有多边形放入四叉树中的一个节点。那么任务就变成了并行,因为每个节点都可以单独处理交集。
这是四叉树的 wiki:
在 3D 中可以使用八叉树。
它实际上什至不必是八叉树。您可以使用任何 space 分区获得相同的结果。你可以找到多边形的最大分离(我们称之为 S
)。并创建说 S/10
space 分区。然后您将有 10 个单独的 space 并行执行。它不仅是并发的,而且不再是 M * N
时间,因为不是每个多边形都必须与其他所有多边形进行比较。
您可以创建多边形与边界框的比率:
这可以在计算上完成一次,以获得平均多边形区域与 BB 比率 R
常数。
或者您可以使用以 BB 为界的圆来处理几何因为您只使用投影多面体:
R = 0.0;
count = 0;
for (each poly) {
count++;
R += polyArea / itsBoundingBoxArea;
}
R = R/count;
然后计算边界框交集的总和。
Sbb = 0.0;
for (box1, box2 where box1.isIntersecting(box2)) {
Sbb += box1.intersect(box2);
}
然后:
Approximation = R * Sbb
如果允许凹多边形,所有这些都将不起作用。因为凹多边形可以占据不到其边界框的 1%。您仍然需要找到凸包。
或者,如果您可以比其外壳更快地找到多边形面积,则可以使用实际计算的平均多边形面积。这也会给你一个不错的近似值,同时避免多边形相交和环绕。