Directamente “trazar” segmentos de línea para numpy array

Uno de mis primeros proyectos realizados en python es la simulación de percolación de palo de Monte Carlo. El código creció continuamente. La primera parte fue la visualización de la percolación del palo. En un área de ancho * largo, una densidad definida (palos / área) de palos rectos con una cierta longitud se traza con coordenadas de inicio aleatorias y dirección. Como a menudo uso gnuplot, escribí las coordenadas de inicio y finalización generadas (x, y) en un archivo de texto para luego hacer gnuplot.

Luego encontré aquí una buena forma de analizar los datos de la imagen usando scipy.ndimage.measurements. La imagen se lee utilizando ndimage.imread en escala de grises. La matriz numpy resultante se reduce aún más a valores booleanos, ya que solo me interesan las conexiones entre diferentes sticks. Los grupos resultantes se analizan con ndimage.measurements. Esto me permite saber si hay vías que se conectan de un lado a otro o no. Un ejemplo minimizado está aquí.

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

Si bien esto funciona principalmente bien, las simulaciones de Monte Carlo con 1000 iteraciones para un mayor número de diferentes densidades de stick terminan funcionando 8 horas o más. Esto se debe en parte al hecho de que las imágenes y matrices creadas son bastante grandes y se trazan miles de palos para densidades más altas. La razón es que quiero simular un rango de geometrías (por ejemplo, longitud entre 500 y 20000 píxeles) mientras minimizo el error debido a la pixelización.

Supongo que la mejor manera sería no usar los datos de imagen y tratarlos como un problema vectorial, pero no tengo ni idea de cómo iniciar un algoritmo. Y las muchas conexiones también pueden dar lugar a grandes matrices de datos.

Manteniendo el método descrito anteriormente, es obvio que escribir datos en un archivo y releerlo no es muy efectivo. Por lo tanto, estoy buscando maneras de acelerar esto. Como primer paso, utilicé matplotlib para crear la imagen, pero al menos al trazar cada stick con una llamada de la ttwig separada, esto es hasta 10 veces más lento para un mayor número de sticks. Crear la lista de coordenadas de barra en una matriz y trazar la lista completa con una sola llamada de plot podría acelerarla pero aún así deja el cuello de botella al escribir y leer la imagen.

¿Puede indicarme un método efectivo para generar directamente la matriz numpy de tipo booleano que representa la imagen en blanco y negro de los palos? ¿Tal vez trazar la lista de coordenadas y convertir la figura de alguna manera a una matriz? También encontré esta interesante discusión donde las líneas se trazan en una imagen PIL. ¿Podría ser esto más rápido que matplotlib?

Dibujar un segmento de línea en una matriz es una capacidad fundamental de cualquier biblioteca de gráficos. El método más simple es probablemente el algoritmo de Bresenham . El algoritmo es simple y rápido: cuando se implementa en un lenguaje rápido, es decir. No recomendaría implementarlo en python puro. Un inconveniente de la versión más simple del algoritmo es que no tiene suavizado. Las líneas muestran “jaggies” . Busque “algoritmos de dibujo de líneas” para obtener métodos más avanzados con mejor suavizado.

Tengo una implementación de Cython del algoritmo de Bresenham en mi paquete eyediagram . La función bres_segment_count incrementa los valores en la matriz de entrada a lo largo de la línea recta desde (x0, y0) hasta (x1, y1). Una modificación que simplemente establece los valores de la matriz en 1 sería un cambio trivial en ese código.

Por ejemplo,

 In [21]: dim = 250 In [22]: num_sticks = 300 

Cada fila de sticks contiene [x0, y0, x1, y1], los puntos finales de un “palo”:

 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 dibuja cada palo usando el algoritmo de Bresenham. Tenga en cuenta que en lugar de simplemente establecer un valor en la línea para, por ejemplo, 1, los valores en img largo de la línea se incrementan.

 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]:  

Aquí está la ttwig que se genera:

trama de palos