【问题标题】：How can I select the pixels that fall within a contour in an image represented by a numpy array?如何在由 numpy 数组表示的图像中选择落在轮廓内的像素？
【发布时间】：2018-11-23 15:52:06
【问题描述】：

我在存储为 2D numpy 数组的图像上绘制了一组轮廓点。轮廓由 x 和 y 坐标的 2 个 numpy 浮点值数组表示。 These coordinates 不是整数，也不能与像素完美对齐，但它们确实可以告诉您轮廓点相对于像素的位置。

我希望能够选择轮廓内的像素。我写了一些与这里给出的答案几乎相同的代码：Access pixel values within a contour boundary using OpenCV in Python

temp_list = []
for a, b in zip(x_pixel_nos, y_pixel_nos):
    temp_list.append([[a, b]]) # 2D array of shape 1x2
temp_array = np.array(temp_list)
contour_array_list = []
contour_array_list.append(temp_array)



lst_intensities = []

# For each list of contour points...
for i in range(len(contour_array_list)):
    # Create a mask image that contains the contour filled in
    cimg = np.zeros_like(pixel_array)
    cv2.drawContours(cimg, contour_array_list, i, color=255, thickness=-1)

# Access the image pixels and create a 1D numpy array then add to list
pts = np.where(cimg == 255)
lst_intensities.append(pixel_array[pts[0], pts[1]])

当我运行这个时，我得到一个错误error: OpenCV(3.4.1) /opt/conda/conda-bld/opencv-suite_1527005509093/work/modules/imgproc/src/drawing.cpp:2515: error: (-215) npoints > 0 in function drawContours

我猜测此时 openCV 对我不起作用，因为我的轮廓是浮点数，而不是整数，openCV 无法处理 drawContours。如果我将轮廓的坐标转换为整数，我会损失很多精度。

那么我怎样才能得到落在轮廓内的像素呢？

这应该是一项微不足道的任务，但到目前为止我还没有找到一种简单的方法来完成它。

【问题讨论】：

如果您的像素是浮点数，那么理论上有无限数量的像素落在轮廓内，因为现在位置是连续的。你想施加什么约束来限制这个像素列表的长度？还要记住，轮廓像素列表是这样的，每个元素都需要是一个具有单一第二维的 3d 数组（即另一篇文章中提到的 N x 1 x 2 数组）。
@rayryeng 我的像素不连续。它们是整数。轮廓是浮动的。我猜这是因为这些轮廓是由一个人在绘制轮廓时放大图像手动绘制的。
啊，我明白了。但是，您的轮廓列表看起来并不正确。请参阅我之前关于它们是 3d 数组的说明。
有趣。我需要进行实验。我现在不在电脑前，所以我要等到今晚晚些时候才能尝试。您是否碰巧在某处有可用的轮廓数据？
@CrisLuengo 目前我并不关心计算效率。如果您能向我展示完成此操作的代码，我将不胜感激。

标签： python numpy image-processing contour opencv-contour

【解决方案1】：

我认为找到所有落在轮廓内的像素的最简单方法如下。

轮廓由一组非整数点描述。我们可以把这些点看成一个多边形的顶点，轮廓就是一个多边形。

我们首先找到多边形的边界框。此边界框之外的任何像素都不在多边形内，不需要考虑。

对于边界框内的像素，我们使用经典测试来测试它们是否在多边形内：从无限远的某个点到该点跟踪一条线，并计算交叉的多边形边（线段）的数量。如果这个数字是奇数，则该点在多边形内。原来Matplotlib包含a very efficient implementation of this algorithm。

我还在习惯 Python 和 Numpy，如果您是 Python 专家，这可能会有点尴尬。但我认为它是直截了当的。首先它计算多边形的边界框，然后创建一个数组points，其中包含该边界框内的所有像素的坐标（我假设像素质心是最重要的）。它将matplotlib.path.contains_points 方法应用于此数组，产生一个布尔数组mask。最后，它重塑这个数组以匹配边界框。

import math
import matplotlib.path
import numpy as np

x_pixel_nos = [...]
y_pixel_nos = [...] # Data from https://gist.github.com/sdoken/173fae1f9d8673ffff5b481b3872a69d

temp_list = []
for a, b in zip(x_pixel_nos, y_pixel_nos):
   temp_list.append([a, b])

polygon = np.array(temp_list)
left = np.min(polygon, axis=0)
right = np.max(polygon, axis=0)
x = np.arange(math.ceil(left[0]), math.floor(right[0])+1)
y = np.arange(math.ceil(left[1]), math.floor(right[1])+1)
xv, yv = np.meshgrid(x, y, indexing='xy')
points = np.hstack((xv.reshape((-1,1)), yv.reshape((-1,1))))

path = matplotlib.path.Path(polygon)
mask = path.contains_points(points)
mask.shape = xv.shape

在此代码之后，需要在图像中定位边界框，并为像素着色。 left 包含图像中与mask 左上角像素对应的像素。

可以提高该算法的性能。如果跟踪来测试像素的光线是水平的，您可以想象沿水平线的所有像素都可以从左侧像素所做的工作中受益。也就是说，计算一条图像线上所有像素的输入/输出状态可能比单个像素的成本稍微多一点。

matplotlib.path.contains_points 算法比对所有点执行单点测试要高效得多，因为适当地对多边形边和顶点进行排序可以使每次测试的成本大大降低，并且在测试多个点时只需进行一次排序立刻。但是这个算法没有考虑到我们要在同一条线上测试很多点。

这些都是我做的时候看到的

pp.plot(x_pixel_nos, y_pixel_nos)
pp.imshow(mask)

使用您的数据运行上述代码后。请注意，y 轴与 imshow 反转，因此是垂直镜像的形状。

【讨论】：

非常感谢您在这里的帮助。我在多个图像上尝试过的这段代码使我能够在“轮廓”中找到像素。（我对代码风格并不挑剔，只关心一个可行的解决方案）我还有几个问题：a）我很想知道在这里使用 hstack 和 meshgrid 时的想法，因为我无法破译你最终是如何使用这些功能的，以及这与它们的正常使用有何不同。 b) 给定一组轮廓点（或顶点），matplotlib.path.Path 是否总能找到凸多边形？如果我们需要的形状是凹的呢？
@itsnotme：多边形不需要是凸的，这段代码中没有这样的假设。 Path 可以编码任何多边形（实际上，我不确定自相交，但你没有这些）。 meshgrid 和 hstack 我以我在 MATLAB 中习惯的方式使用。基本上，meshgrid 在xv 和yv 中生成二维数组，它们一起包含矩形区域内每个点的 (x,y) 坐标，由向量 x 和 y 描述。正如contains_points 所期望的那样，我将这两个数组重新整形为一维数组并将它们连接起来（hstack）以形成一个 2xN 坐标数组。
在这里，我绘制了来自您的代码的轮廓和像素。 github.com/sdoken/pydicom_rtstruct_contour_mapping/blob/master/…如果你仔细观察彩色区域的底部，你会发现绿色区域缺少底部图像中红色边界的凹度。我不确定这是一个实际问题（算法选择了一些在没有的点时处于其中的点）还是可视化过程的产物，而不是问题
@itsnotme：请参阅我刚刚附在答案中的图表。没有理由相信被称为创建凸包的函数。我不确定您发布的图片中发生了什么。该多边形是否有任何自相交？

【解决方案2】：

借助 python 中的 Shapely 库，可以轻松做到：

从 shapely.geometry 导入点、多边形

将所有 x,y 坐标转换为形状优美的多边形：

coords = [(0, 0), (0, 2), (1, 1), (2, 2), (2, 0), (1, 1), (0, 0)]
pl = Polygon(coords)

现在在每个多边形中查找像素为：

minx, miny, maxx, maxy = pl.bounds
minx, miny, maxx, maxy = int(minx), int(miny), int(maxx), int(maxy)
box_patch = [[x,y] for x in range(minx,maxx+1) for y in range(miny,maxy+1)]
pixels = []
for pb in box_patch: 
  pt = Point(pb[0],pb[1])
  if(pl.contains(pt)):
    pixels.append([int(pb[0]), int(pb[1])])
return pixels

将此循环用于每组坐标，然后用于每个多边形。

很高兴 :)

【讨论】：

【解决方案3】：

skimage.draw.polygon 可以处理此1，请参阅该页面上此函数的示例代码。

如果你只想要轮廓，你可以做skimage.segmentation.find_boundaries2。

【讨论】：