有一个整数矩阵 MxN 如何将它们分组为具有增强几何形状的多边形？
Posted
技术标签:
【中文标题】有一个整数矩阵 MxN 如何将它们分组为具有增强几何形状的多边形？【英文标题】：Having a matrix MxN of integers how to group them into polygons with boost geometry? 【发布时间】：2011-11-07 17:05:52 【问题描述】：我们有一个给定整数的矩阵（从 1 到 INT_MAX 的任意一个）
1 2 3
1 3 3
1 3 3
100 2 1
我们希望为矩阵中的每个唯一 int 创建具有相同颜色的多边形，因此我们的多边形将具有如下所示的坐标/分组。
我们可以生成这样的图像：
哪个*（因为执行的vectirisation会缩放到这样的大小）：
（很糟糕的图见谅）
是否有可能以及如何使用 boost 几何做这样的事情？
更新：
所以@sehe sad: I'd simply let Boost Geometry do most of the work. 所以我使用纯 Boost.Geometry 逐个像素地创建了 this 像素类 aeria 种植者，编译、运行但我需要它在集群数据上运行.. 我有 1000 x 1800 个 uchars 文件（每个唯一的 uchar == 数据属于该分类器）。这段代码的问题：在第 18 行，它变得非常缓慢，以至于每个点创建开始需要超过一秒=(
代码：
//Boost
#include <boost/assign.hpp>
#include <boost/foreach.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/algorithm/string.hpp>
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/geometries.hpp>
#include <boost/geometry/multi/geometries/multi_polygon.hpp>
#include <boost/geometry/geometries/adapted/boost_tuple.hpp>
//and this is why we use Boost Geometry from Boost trunk
//#include <boost/geometry/extensions/io/svg/svg_mapper.hpp>
BOOST_GEOMETRY_REGISTER_BOOST_TUPLE_CS(cs::cartesian)
void make_point(int x, int y, boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > & ring)
using namespace boost::assign;
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y-1));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x, y-1));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x, y));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y-1));
boost::geometry::correct(ring);
void create_point(int x, int y, boost::geometry::model::multi_polygon< boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > > & mp)
boost::geometry::model::multi_polygon< boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > > temp;
boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > ring;
make_point(x, y, ring);
boost::geometry::union_(mp, ring, temp);
boost::geometry::correct(temp);
mp=temp;
int main()
using namespace boost::assign;
boost::geometry::model::multi_polygon< boost::geometry::model::polygon < boost::geometry::model::d2::point_xy<double> > > pol, simpl;
//read image
std::ifstream in("1.mask", std::ios_base::in | std::ios_base::binary);
int sx, sy;
in.read(reinterpret_cast<char*>(&sy), sizeof(int));
in.read(reinterpret_cast<char*>(&sx), sizeof(int));
std::vector< std::vector<unsigned char> > image(sy);
for(int i =1; i <= sy; i++)
std::vector<unsigned char> row(sx);
in.read(reinterpret_cast<char*>(&row[0]), sx);
image[i-1] = row;
//
std::map<unsigned char, boost::geometry::model::multi_polygon < boost::geometry::model::polygon < boost::geometry::model::d2::point_xy<double> > > > layered_image;
for(int y=1; y <= sy; y++)
for(int x=1; x <= sx; x++)
if (image[y-1][x-1] != 1)
create_point(x, y, layered_image[image[y-1][x-1]]);
std::cout << x << " : " << y << std::endl;
所以你可以看到我的代码 suks.. 所以我决定为 @sehe code 创建一个渲染器：
#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include <set>
//Boost
#include <boost/assign.hpp>
#include <boost/array.hpp>
#include <boost/algorithm/string.hpp>
#include <boost/geometry/geometry.hpp>
#include <boost/geometry/geometries/geometries.hpp>
#include <boost/geometry/multi/geometries/multi_polygon.hpp>
#include <boost/geometry/geometries/adapted/boost_tuple.hpp>
#include <boost/random/uniform_int.hpp>
#include <boost/random/mersenne_twister.hpp>
//and this is why we use Boost Geometry from Boost trunk
#include <boost/geometry/extensions/io/svg/svg_mapper.hpp>
BOOST_GEOMETRY_REGISTER_BOOST_TUPLE_CS(cs::cartesian)
namespace mxdetail
typedef size_t cell_id; // row * COLS + col
template <typename T> struct area
T value;
typedef std::vector<cell_id> cells_t;
cells_t cells;
;
template <typename T, size_t Rows, size_t Cols>
std::vector<area<T> > getareas(const boost::array<boost::array<T, Cols>, Rows>& matrix)
typedef boost::array<boost::array<T, Cols>, Rows> mtx;
std::vector<area<T> > areas;
struct visitor_t
const mtx& matrix;
std::set<cell_id> visited;
visitor_t(const mtx& mtx) : matrix(mtx)
area<T> start(const int row, const int col)
area<T> result;
visit(row, col, result);
return result;
void visit(const int row, const int col, area<T>& current)
const cell_id id = row*Cols+col;
if (visited.end() != visited.find(id))
return;
bool matches = current.cells.empty() || (matrix[row][col] == current.value);
if (matches)
visited.insert(id);
current.value = matrix[row][col];
current.cells.push_back(id);
// process neighbours
for (int nrow=std::max(0, row-1); nrow < std::min((int) Rows, row+2); nrow++)
for (int ncol=std::max(0, col-1); ncol < std::min((int) Cols, col+2); ncol++)
/* if (ncol!=col || nrow!=row) */
visit(nrow, ncol, current);
visitor(matrix);
for (int r=0; r < (int) Rows; r++)
for (int c=0; c < (int) Cols; c++)
mxdetail::area<int> area = visitor.start(r,c);
if (!area.cells.empty()) // happens when startpoint already visited
areas.push_back(area);
return areas;
typedef boost::array<int, 4> row;
template <typename T, size_t N>
boost::array<T, N> make_array(const T (&a)[N])
boost::array<T, N> result;
std::copy(a, a+N, result.begin());
return result;
void make_point(int x, int y, boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > & ring)
using namespace boost::assign;
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y-1));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x, y-1));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x, y));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y));
boost::geometry::append( ring, boost::geometry::model::d2::point_xy<double>(x-1, y-1));
boost::geometry::correct(ring);
void create_point(int x, int y, boost::geometry::model::multi_polygon< boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > > & mp)
boost::geometry::model::multi_polygon< boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > > temp;
boost::geometry::model::polygon<boost::geometry::model::d2::point_xy<double> > ring;
make_point(x, y, ring);
boost::geometry::union_(mp, ring, temp);
boost::geometry::correct(temp);
mp=temp;
boost::random::mt19937 rng;
boost::random::uniform_int_distribution<> color(10,255);
std::string fill_rule()
int red, green, blue;
red = color(rng);
green = color(rng);
blue = color(rng);
std::ostringstream rule;
rule << "fill-rule:nonzero;fill-opacity:0.5;fill:rgb("
<< red << "," << green << "," << blue
<< ");stroke:rgb("
<< (red - 5) << "," << (green - 5) << "," << (blue -5)
<< ");stroke-width:2";
return rule.str();
int main()
int sx = 4;
int sy = 5;
int row0[] = 1 , 2, 3, 3, ;
int row1[] = 1 , 3, 3, 3,;
int row2[] = 1 , 3, 3, 3, ;
int row3[] = 2 , 2, 1, 2, ;
int row4[] = 100, 2, 2, 2, ;
boost::array<row, 5> matrix;
matrix[0] = make_array(row0);
matrix[1] = make_array(row1);
matrix[2] = make_array(row2);
matrix[3] = make_array(row3);
matrix[4] = make_array(row4);
typedef std::vector<mxdetail::area<int> > areas_t;
typedef areas_t::value_type::cells_t cells_t;
areas_t areas = mxdetail::getareas(matrix);
using namespace boost::assign;
typedef boost::geometry::model::polygon
<
boost::geometry::model::d2::point_xy<double>
> polygon;
typedef boost::geometry::model::multi_polygon<polygon> mp;
typedef boost::geometry::point_type<mp>::type point_type;
std::string filename = "draw.svg";
std::ofstream svg(filename.c_str());
boost::geometry::svg_mapper<point_type> mapper(svg, 400, 400);
for (areas_t::const_iterator it=areas.begin(); it!=areas.end(); ++it)
mp pol;
std::cout << "area of " << it->value << ": ";
for (cells_t::const_iterator pit=it->cells.begin(); pit!=it->cells.end(); ++pit)
int row = *pit / 3, col = *pit % 3;
std::cout << "(" << row << "," << col << "), ";
create_point( (row+1), (col+1), pol);
std::cout << std::endl;
mapper.add(pol);
mapper.map(pol, fill_rule());
std::cout << "areas detected: " << areas.size() << std::endl;
std::cin.get();
这段代码是可编译的，但它很烂（看来我毕竟不知道如何使用数组......）：
【问题讨论】：
这很有趣。您能否详细说明“我们的多边形将具有坐标/分组”的含义。在您的示例中，坐标是什么？你的意思是，行/列索引就是坐标？ 相关：meta.stackexchange.com/questions/19478/the-many-memes-of-meta/… 【参考方案1】：简而言之，如果我的问题正确，我会让 Boost Geometry 完成大部分工作。
对于 NxM 的样本矩阵，创建 NxM 'flyweight' 矩形多边形以对应于每个矩阵单元视觉上。
现在，使用迭代加深算法，找到所有组：
* for each _unvisited_ cell in matrix
* start a new group
* [visit:]
- mark _visited_
- for each neighbour with equal value:
- add to curent group and
- recurse [visit:]
请注意，此算法的结果可能是具有相同值的不同组（表示不相交的多边形）。例如。 OP 中样本中的值 2 将导致 两个 组。
现在对于每个组，只需调用 Boost Geometry's Union_ algorithm 即可找到表示该组的合并多边形。
示例实现
更新这是 C++11 中未优化的实现：
编辑 查看这里C++03 version (using Boost)
测试中使用的样本数据对应于问题中的矩阵。
#include <iostream>
#include <array>
#include <vector>
#include <set>
namespace mxdetail
typedef size_t cell_id; // row * COLS + col
template <typename T> struct area
T value;
std::vector<cell_id> cells;
;
template <typename T, size_t Rows, size_t Cols>
std::vector<area<T> > getareas(const std::array<std::array<T, Cols>, Rows>& matrix)
typedef std::array<std::array<T, Cols>, Rows> mtx;
std::vector<area<T> > areas;
struct visitor_t
const mtx& matrix;
std::set<cell_id> visited;
visitor_t(const mtx& mtx) : matrix(mtx)
area<T> start(const int row, const int col)
area<T> result;
visit(row, col, result);
return result;
void visit(const int row, const int col, area<T>& current)
const cell_id id = row*Cols+col;
if (visited.end() != visited.find(id))
return;
bool matches = current.cells.empty() || (matrix[row][col] == current.value);
if (matches)
visited.insert(id);
current.value = matrix[row][col];
current.cells.push_back(id);
// process neighbours
for (int nrow=std::max(0, row-1); nrow < std::min((int) Rows, row+2); nrow++)
for (int ncol=std::max(0, col-1); ncol < std::min((int) Cols, col+2); ncol++)
/* if (ncol!=col || nrow!=row) */
visit(nrow, ncol, current);
visitor(matrix);
for (int r=0; r < Rows; r++)
for (int c=0; c < Cols; c++)
auto area = visitor.start(r,c);
if (!area.cells.empty()) // happens when startpoint already visited
areas.push_back(area);
return areas;
int main()
typedef std::array<int, 3> row;
std::array<row, 4> matrix =
row 1 , 2, 3, ,
row 1 , 3, 3, ,
row 1 , 3, 3, ,
row 100, 2, 1, ,
;
auto areas = mxdetail::getareas(matrix);
std::cout << "areas detected: " << areas.size() << std::endl;
for (const auto& area : areas)
std::cout << "area of " << area.value << ": ";
for (auto pt : area.cells)
int row = pt / 3, col = pt % 3;
std::cout << "(" << row << "," << col << "), ";
std::cout << std::endl;
使用gcc-4.6 -std=c++0x编译输出为：
areas detected: 6
area of 1: (0,0), (1,0), (2,0),
area of 2: (0,1),
area of 3: (0,2), (1,1), (1,2), (2,1), (2,2),
area of 100: (3,0),
area of 2: (3,1),
area of 1: (3,2),
【讨论】：
所以我开始学习它......并遇到了一个小问题。描述它here (with compilable source code) 添加了我的一些代码...请注意 - 纯提升几何...)) 以及对问题的更深入描述...终于加入帐户 - Kiss Kabumbus Blender(aka Ole Jak).. . 发现每月有 30 个 Q 限制=( Got 5k rep...=) @Kabumbus：好吧，那么也许是你发布了这个问题？我刚刚发布了我的示例代码的仅限 C++03 版本（并修复了一个错误）。见***.com/questions/8067038/… 已创建渲染器...尝试稍微扩展原始数组，但在生成的图像中遇到了一些奇怪的问题...请参阅更新后的所有代码 有趣。我会看看我是否可以给予它一些关注。但是，这可能会在明天晚上 - 很晚。【参考方案2】：当点数很大时（例如，超过 1000x1000），上面的解决方案会占用大量内存。这正是 topic-starter 发生的事情。
下面我展示了更具可扩展性的方法。
这里我将两个问题分开：一个是找到区域，另一个是将它们转换为多边形。
第一个问题实际上等同于找到网格图的连通分量，其中邻居有边当且仅当它们具有相同的“颜色”附加到它上面。可以使用来自boost-graph 的网格图。
#include <boost/graph/grid_graph.hpp>
// Define dimension lengths, a MxN in this case
boost::array<int, 2> lengths = M, N ;
// Create a MxN two-dimensional, unwrapped grid graph
grid_graph<2> graph(lengths);
接下来，我们应该将给定的矩阵 M 转换为边缘过滤器：如果邻居的“颜色”相同，则存在网格边缘。
template <class Matrix>
struct GridEdgeFilter
typedef grid_graph<2> grid;
GridEdgeFilter(const Matrix & m, const grid&g):_m(&m),_g(&g)
/// \return true iff edge is present in the graph
bool operator()(grid::edge_descriptor e) const
grid::vertex_descriptor src = source(e,*_g), tgt = target(e,*_g);
//src[0] is x-coord of src, etc. The value (*m)[x,y] is the color of the point (x,y).
//Edge is preserved iff matrix values are equal
return (*_m)[src[0],src[1]] == (*_m)[tgt[0],tgt[1]];
const Matrix * _m;
const grid* _g;
;
最后，我们定义了一个boost::filtered_graph的grid和EdgeFilter，并为连通分量调用Boost.Graph算法。
每个连接的组件代表一组单一颜色的点，即我们想要转换为多边形的区域。
这里我们有另一个问题。 Boost.Geometry 只允许一个一个地合并多边形。因此，当多边形数量很大时，它会变得非常慢。
更好的方法是使用 Boost.Polygon，即它的Property Merge 功能。一个从空的property_merge 对象开始，然后插入给定颜色的矩形（您可以将颜色设置为属性）。然后调用方法merge 并获得polygon_set 作为输出。
【讨论】：
以上是关于有一个整数矩阵 MxN 如何将它们分组为具有增强几何形状的多边形？的主要内容，如果未能解决你的问题，请参考以下文章