直接"plot"线段转成numpy数组

Question

我在 python 中实现的第一个项目之一是 Monte Carlo 棒渗透模拟。代码不断增长。第一部分是棍子渗滤的可视化。在宽度*长度的区域中，使用随机起始坐标和方向绘制具有一定长度的直杆的定义密度（sticks/area）。因为我经常使用 gnuplot，所以我将生成的 (x, y) 开始和结束坐标写入一个文本文件，以便之后进行 gnuplot。

然后我发现 here 使用 scipy.ndimage.measurements 分析图像数据的好方法。使用灰度 ndimage.imread 读取图像。生成的 numpy 数组进一步简化为布尔值，因为我只对不同棒之间的连接感兴趣。然后使用 ndimage.measurements 分析生成的聚类。这使我能够确定是否存在从一侧连接到另一侧的路径。此处为最小化示例。

import random
import math
from scipy.ndimage import measurements
from scipy.ndimage import imread
import numpy as np
import matplotlib.pyplot as plt

#dimensions of plot
width = 10
length = 8
stick_length = 1
fig = plt.figure(frameon=False)
ax = fig.add_axes([0, 0, 1, 1])
fig.set_figwidth(width)
fig.set_figheight(length)
ax.axis('off')

file = open("coordinates.txt", "w")

for i in range (300):
    # randomly create (x,y) start coordinates in channel and direction
    xstart = width * random.random()        # xstart = 18
    ystart = length * random.random()        # ystart = 2
    # randomly generate direction of stick from start coordinates and convert from GRAD in RAD
    dirgrad = 360 * random.random()
    dirrad = math.radians(dirgrad)
    # calculate (x,y) end coordinates
    xend = xstart + (math.cos(dirrad) * stick_length)
    yend = ystart + (math.sin(dirrad) * stick_length)
    # write start and end coordinates into text file for gnuplot plotting
    file.write(str(i) + ":\t" + str(xstart) + "\t" + str(ystart) + "\t" + str(dirgrad) + ":\t" + str(xend) + "\t" + str(yend) + "\n")
    file.write(str(i) + ":\t" + str(xend) + "\t" + str(yend) + "\n\n")
    # or plot directly with matplotlib
    ax.plot([xstart,xend],[ystart,yend],"black", lw=1)
fig.savefig("testimage.png", dpi=100)

# now read just saved image and do analysis with scipy.ndimage
fig1, ax1 = plt.subplots(1,1)
img_input = imread("testimage.png", flatten = True) # read image to np.ndarray in grey scales
img_bw = img_input < 255 # convert grey scales to b/w (boolean)
labeled_array, num_clusters = measurements.label(img_bw) #labeled_array: labeled clusters in array, num_clusters: number of clusters
area = measurements.sum(img_bw, labeled_array, index=np.arange(labeled_array.max() + 1)) # area of each cluster
areaImg = area[labeled_array] # label each cluster with labelnumber=area
cax = ax1.imshow(areaImg, origin='upper', interpolation='nearest', cmap = 'rainbow')
cbar = fig1.colorbar(cax)
fig1.savefig("testimage_analyzed.png")

虽然这主要工作得很好 Monte Carlo 对大量不同棒密度进行 1000 次迭代的模拟最终运行 8 小时或更长时间。这部分是由于这样一个事实，即创建的图像和阵列非常大，并且为了更高的密度绘制了数千根棒。原因是我想模拟一系列几何形状（例如长度在 500 到 20000 像素之间），同时最大限度地减少由于像素化引起的误差。

我想最好的方法是不使用图像数据并将其视为矢量问题，但我什至不知道如何开始算法。许多连接也可能导致大型数据阵列。

继续使用上述方法，很明显将数据写入文件并重新读取它并不是很有效。因此，我正在寻找加快速度的方法。作为第一步，我使用 matplotlib 来创建图像，但是至少在使用单独的绘图调用绘制每根棍子时，对于更多的棍子，这会慢 10 倍。在数组中创建棒坐标列表并通过一次绘图调用绘制完整列表可能会加快速度，但仍然存在写入和读取图像的瓶颈。

你能指点一个有效的方法来直接生成代表木棍黑白图像的布尔类型numpy数组吗？也许绘制坐标列表并以某种方式将图形转换为数组？我还发现这个有趣的 discussion，其中线条被绘制到 PIL 图像上。这可能比 matplotlib 更快吗？

Answer 1

在数组中绘制线段是任何图形库的基本功能。最简单的方法大概是Bresenham's algorithm. The algorithm is simple and fast--when implemented in a fast language, that is. I wouldn't recommend implementing it in pure python. A drawback of the simplest version of the algorithm is that it is not anti-aliased. The lines show "jaggies"。搜索 "line drawing algorithms" 以获得具有更好抗锯齿效果的更高级方法。

我有一个Cython implementation of Bresenham's algorithm in my eyediagram package。函数 bres_segment_count 沿着从 (x0, y0) 到 (x1, y1) 的直线递增输入数组中的值。仅将数组值设置为 1 的修改对该代码的改动很小。

例如，

In [21]: dim = 250

In [22]: num_sticks = 300

sticks的每一行包含[x0,y0,x1,y1]，"stick"的端点：

In [23]: sticks = np.random.randint(0, dim, size=(num_sticks, 4)).astype(np.int32)

In [24]: img = np.zeros((dim, dim), dtype=np.int32)

bres_segments_count 使用 Bresenham 算法绘制每根木棍。请注意，不是简单地将行中的值设置为 1，而是递增 img 中沿行的值。

In [25]: from eyediagram._brescount import bres_segments_count

In [26]: bres_segments_count(sticks, img)

In [27]: plt.imshow(img, interpolation='nearest', cmap=cm.hot)
Out[27]: <matplotlib.image.AxesImage at 0x10f94b110>

这是生成的图：

直接"plot"线段转成numpy数组

Directly "plot" line segments to numpy array

python

matplotlib

scipy

montecarlo

ndimage