我在python中實現的第一個項目是Monte Carlo模擬棒狀滲透。代碼不斷增長。第一部分是棒滲透的可視化。在寬度*長度的區域中,用一定長度的直杆的定義密度(棒/面積)用隨機開始座標和方向繪製。由於我經常使用gnuplot,因此我將生成的(x,y)開始和結束座標寫入文本文件,以便在之後進行gnuplot繪圖。直接「繪製」線段到numpy陣列
然後,我發現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快嗎?
謝謝你指點我到這個方向。在我首次在純python中實現它之前,沒有使用過cython,這已經有了改進。然而,我現在正在努力讓你的代碼運行。下載完整的軟件包並保存在「cython」目錄中。運行時 – MrCyclophil
最後在Linux下運行。不知道從哪裏獲得Visual Studio 2010 Express以使其在Windows下運行。 – MrCyclophil
僅供參考,如果有人在Win8.1 64bit上安裝C編譯器時也有問題。對於Python 3.5,它變得更容易。只需下載並安裝免費的Visual Studio 2015社區版。 [本博客](http://blog.ionelmc.ro/2014/12/21/compiling-python-extensions-on-windows/)爲Python 2.7,3.4和3.5提供了一個很好的總結。 – MrCyclophil