1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
|
import numpy as np
import pycuda.driver as cuda
from gputhread import GPUThread
from detectors import minilbne
from Queue import Queue
from threadtest import create_job
import matplotlib.pyplot as plt
from itertoolset import roundrobin
from color import map_wavelength
from solids import r7081
from geometry import Geometry
nphotons = 1000
jobs = Queue()
output = Queue()
geometry = minilbne
geometry.build(bits=8)
cuda.init()
try:
gputhread = GPUThread(5, geometry, jobs, output, 64)
gputhread.start()
job = create_job((0,0,0), nphotons)
x = np.empty((nphotons, 10, 3))
for i in range(10):
print '%i' % i
x[:,i,0] = job.positions['x']
x[:,i,1] = job.positions['y']
x[:,i,2] = job.positions['z']
jobs.put(job)
jobs.join()
job = output.get()
for j in np.unique(job.states):
print 'state %2i, %i' % (j, len(job.states[job.states == j]))
finally:
gputhread.stop()
gputhread.join()
mask = job.states == 2
rgb = (map_wavelength(job.wavelengths[mask])*255).astype(np.uint32)
def format_hex_string(s):
return '#' + s.rstrip('L')[2:].zfill(6)
colors = map(format_hex_string, map(hex, rgb[:,0] << 16 | rgb[:,1] << 8 | rgb[:,2]))
plt.figure()
plt.plot(*roundrobin(x[mask,:,0], x[mask,:,2], colors))
plt.show()
|
