直接“绘制"图2.线段到numpy数组 [英] Directly "plot" line segments to numpy array

问题描述

我在python中实现的第一个项目之一是对棒渗透进行了Monte Carlo模拟.代码不断增长. 第一部分是棒渗滤的可视化.在宽度×长度的区域中，使用随机的起始坐标和方向绘制具有一定长度的直棒的定义密度(棒/面积).当我经常使用gnuplot时，我将生成的(x，y)起始和结束坐标写入文本文件中，然后对其进行gnuplot.

One of my first projects realized in python does Monte Carlo simulation of stick percolation. The code grew continually. The first part was the visualization of the stick percolation. In an area of width*length a defined density (sticks/area) of straight sticks with a certain length are plotted with random start coordinates and direction. As I often use gnuplot I wrote the generated (x, y) start and end coordinates to a text file to gnuplot them afterwards.

然后，我找到了此处一种使用scipy.ndimage.measurements分析图像数据的好方法.通过使用ndimage.imread灰度读取图像.所得的numpy数组被进一步简化为布尔值，因为我只对不同棍子之间的连接感兴趣.然后使用ndimage.measurements分析所得的群集.这使我可以发现是否存在从一侧连接到另一侧的路径.一个最小的示例在这里.

I then found here a nice way to analyze the image data using scipy.ndimage.measurements. The image is read by using ndimage.imread in greyscales. The resulting numpy array is further reduced to boolean values as I am only interested in connections between different sticks. The resulting clusters are then analyzed with ndimage.measurements. This allows me to find out if there are pathways that connect from one side to the other or not. A minimized example is here.

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")

虽然这基本上可以正常工作，但对于大量不同的杆密度，经过1000次迭代的精细蒙特卡洛模拟最终运行了8个小时或更长时间.这部分是由于这样的事实，即所创建的图像和阵列非常大，并且绘制了成千上万的棒以实现更高的密度.原因是我想模拟一系列几何图形(例如，长度在500到20000像素之间)，同时最大程度地减少由于像素化引起的误差.

While this works principally just fine Monte Carlo simulations with 1000 iterations for a larger number of different stick densities end up running 8 hours or more. This is partly due to the fact, that the created images & arrays are quite large and thousands of sticks are plotted for higher densities. The reason is that I want to simulate a range of geometries (eg length between 500 and 20000 pixels) while minimizing the error due to the pixelization.

我猜最好的方法是不使用图像数据并将其视为向量问题，但是我什至不知道如何启动算法.而且许多连接也可能导致大型数据数组.

I guess the best way would be not to use image data and treat it as a vector problem, but I have no idea how to even start an algorithm. And the many connections might result in large data arrays as well.

使用上述方法，很明显将数据写入文件并重新读取不是很有效.因此，我正在寻找加快速度的方法.第一步，我使用matplotlib来创建图像，但是至少在使用单独的绘图调用来绘制每个棍子时，对于大量的棍子，这要慢10倍.在数组中创建摇杆坐标列表并使用一次绘图调用绘制完整列表可能会加快速度，但仍然留下了编写和读取图像的瓶颈.

Staying with the above described method it is obvious that writing data to a file and re-reading it is not very effective. I am therefore looking for ways to speed this up. As a first step I used matplotlib to create the image however at least when plotting each stick with a separate plot call this is up to 10x slower for a larger number of sticks. Creating the list of stick coordinates in an array and plotting the complete list with one plot call might speed it up but still leaves the bottleneck of writing and reading the image.

您能指出一个有效的方法来直接生成表示棒的黑白图像的布尔型numpy数组吗?也许绘制坐标列表并将图形以某种方式转换为数组?我还发现了有趣的讨论，其中绘制了线条到PIL图像.会比matplotlib快吗?

Can you point me to an effective method to directly generate the boolean type numpy array representing the black and white image of the sticks? Maybe plot the list of coordinates and convert the figure in some way to an array? I also found this interesting discussion where lines are plotted to a PIL image. Might this be faster than matplotlib?

推荐答案

在数组中绘制线段是任何图形库的基本功能.最简单的方法可能是 Bresenham算法.该算法既简单又快速-即以快速语言实现时.我不建议在纯python中实现它.该算法最简单版本的一个缺点是它没有抗锯齿.这些行显示"jaggies" .搜索线条绘制算法"以获取具有更好抗锯齿功能的更高级方法.

Drawing a line segment in an array is a fundamental capability of any graphics library. The simplest method is probably 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". Search for "line drawing algorithms" for more advanced methods with better anti-aliasing.

我有一个 Bresenham算法的Cython实现我的眼图包.函数bres_segment_count沿直线将输入数组中的值从(x0，y0)递增到(x1，y1).仅将数组值设置为1的修改将是对该代码的微不足道的更改.

I have a Cython implementation of Bresenham's algorithm in my eyediagram package. The function bres_segment_count increments the values in the input array along the straight line from (x0, y0) to (x1, y1). A modification that simply sets the array values to 1 would be a trivial change to that code.

例如，

In [21]: dim = 250

In [22]: num_sticks = 300

sticks的每一行都包含[x0，y0，x1，y1]，即棒"的端点:

Each row of sticks contains [x0, y0, x1, y1], the end points of a "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使用布雷森纳姆算法绘制每个棍子.请注意，不是简单地将行中的值设置为1，而是沿行中的img中的值增加.

bres_segments_count draws each stick using Bresenham's algorithm. Note that instead of simply setting a value in the line to, say, 1, the values in img along the line are incremented.

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>

这是生成的图:

这篇关于直接“绘制"图2.线段到numpy数组的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！