summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--chroma/camera.py2
-rw-r--r--chroma/cuda/physical_constants.h6
-rw-r--r--chroma/demo/__init__.py13
-rw-r--r--chroma/demo/checkerboard.py19
-rw-r--r--chroma/demo/optics.py6
-rw-r--r--chroma/demo/pmt.py2
-rw-r--r--chroma/generator/g4gen.py4
-rw-r--r--chroma/gpu/daq.py2
-rw-r--r--chroma/gpu/tools.py1
-rw-r--r--chroma/likelihood.py4
-rw-r--r--chroma/pmt.py7
-rw-r--r--chroma/project.py2
-rw-r--r--src/G4chroma.cc4
-rw-r--r--test/data/ray_intersection.npybin0 -> 1920080 bytes
-rw-r--r--test/data/ray_intersection.npzbin1920080 -> 0 bytes
-rw-r--r--test/test_ray_intersection.py6
16 files changed, 39 insertions, 39 deletions
diff --git a/chroma/camera.py b/chroma/camera.py
index 32c7e88..628307a 100644
--- a/chroma/camera.py
+++ b/chroma/camera.py
@@ -143,7 +143,7 @@ class Camera(multiprocessing.Process):
self.axis1 = np.array([0,0,1], float)
self.axis2 = np.array([1,0,0], float)
- self.film_width = 0.035
+ self.film_width = 35.0 # mm
pos, dir = from_film(self.point, axis1=self.axis1, axis2=self.axis2,
size=self.size, width=self.film_width)
diff --git a/chroma/cuda/physical_constants.h b/chroma/cuda/physical_constants.h
index 2ff87cd..0082503 100644
--- a/chroma/cuda/physical_constants.h
+++ b/chroma/cuda/physical_constants.h
@@ -1,7 +1,9 @@
#ifndef __PHYSICAL_CONSTANTS_H__
#define __PHYSICAL_CONSTANTS_H__
-#define SPEED_OF_LIGHT 2.99792458e8
-#define PI 3.141592653589793
+// mm/ns
+#define SPEED_OF_LIGHT 299.792458f
+
+#define PI 3.141592653589793f
#endif
diff --git a/chroma/demo/__init__.py b/chroma/demo/__init__.py
index edd344a..7df7df7 100644
--- a/chroma/demo/__init__.py
+++ b/chroma/demo/__init__.py
@@ -7,10 +7,10 @@ import numpy.linalg
from chroma.make import sphere
from chroma.stl import mesh_from_stl
from chroma.geometry import *
-from chroma.transform import rotate, make_rotation_matrix
+from chroma.transform import rotate, make_rotation_matrix, normalize
from chroma.demo.pmt import build_8inch_pmt_with_lc
-from chroma.demo.optics import water
+from chroma.demo.optics import water, black_surface
from chroma.demo.checkerboard import build_checkerboard_scene as checkerboard_scene
@@ -24,22 +24,23 @@ def spherical_spiral(radius, spacing):
while t < np.pi:
yield np.array([sin(t) * sin(a*t), sin(t) * cos(a*t), cos(t)])*radius
- dt = spacing / sqrt(1 + a**2 * sin(t) ** 2)
+ dt = dl / sqrt(1 + a**2 * sin(t) ** 2)
t += dt
-def detector(pmt_radius=30.0, sphere_radius=30.5, spiral_step=0.14):
+def detector(pmt_radius=14000.0, sphere_radius=14500.0, spiral_step=350.0):
pmt = build_8inch_pmt_with_lc()
geo = Geometry(water)
geo.add_solid(Solid(sphere(sphere_radius,nsteps=200),
water, water,
+ surface=black_surface,
color=0xBBFFFFFF))
geo.pmtids = []
for position in spherical_spiral(pmt_radius, spiral_step):
- direction = -position/numpy.linalg.norm(position)
+ direction = -normalize(position)
# Orient PMT that starts facing Y axis
y_axis = np.array((0.0,1.0,0.0))
@@ -55,5 +56,5 @@ def detector(pmt_radius=30.0, sphere_radius=30.5, spiral_step=0.14):
return geo
def tiny():
- return detector(1.0, 1.5, 0.4)
+ return detector(2000.0, 2500.0, 700.0)
diff --git a/chroma/demo/checkerboard.py b/chroma/demo/checkerboard.py
index dd2e3e1..801f9da 100644
--- a/chroma/demo/checkerboard.py
+++ b/chroma/demo/checkerboard.py
@@ -6,8 +6,8 @@ from chroma.make import sphere
from chroma.demo.optics import *
def build_checkerboard_scene(checkers_per_side=10, squares_per_checker=50):
- x = np.linspace(-5.0, 5.0, checkers_per_side*squares_per_checker+1)
- y = np.linspace(-5.0, 5.0, checkers_per_side*squares_per_checker+1)
+ x = np.linspace(-5000.0, 5000.0, checkers_per_side*squares_per_checker+1)
+ y = np.linspace(-5000.0, 5000.0, checkers_per_side*squares_per_checker+1)
vertices = np.array(tuple(product(x,y,[0])))
@@ -29,14 +29,15 @@ def build_checkerboard_scene(checkers_per_side=10, squares_per_checker=50):
checkerboard = Solid(checkerboard_mesh, vacuum, vacuum, surface=checkerboard_surface, color=checkerboard_color)
- sphere1 = Solid(sphere(nsteps=512), water, vacuum)
- sphere2 = Solid(sphere(nsteps=512), vacuum, vacuum, surface=shiny_surface)
- sphere3 = Solid(sphere(nsteps=512), vacuum, vacuum, surface=lambertian_surface)
+ sphere1 = Solid(sphere(1000.0, nsteps=512), water, vacuum)
+ sphere2 = Solid(sphere(1000.0, nsteps=512), vacuum, vacuum,
+ surface=shiny_surface)
+ sphere3 = Solid(sphere(1000.0, nsteps=512), vacuum, vacuum, surface=lambertian_surface)
checkerboard_scene = Geometry()
- checkerboard_scene.add_solid(checkerboard, displacement=(0,0,-1.5))
- checkerboard_scene.add_solid(sphere1, displacement=(2.0,-2.0,0))
- checkerboard_scene.add_solid(sphere2, displacement=(-2.0,-2.0,0))
- checkerboard_scene.add_solid(sphere3, displacement=(0.0,2.0,0))
+ checkerboard_scene.add_solid(checkerboard, displacement=(0,0,-1500.0))
+ checkerboard_scene.add_solid(sphere1, displacement=(2000.0,-2000.0,0))
+ checkerboard_scene.add_solid(sphere2, displacement=(-2000.0,-2000.0,0))
+ checkerboard_scene.add_solid(sphere3, displacement=(0.0,2000.0,0))
return checkerboard_scene
diff --git a/chroma/demo/optics.py b/chroma/demo/optics.py
index 2fd34d0..db14032 100644
--- a/chroma/demo/optics.py
+++ b/chroma/demo/optics.py
@@ -51,7 +51,7 @@ r7081hqe_photocathode.set('reflect_diffuse', 1.0 - r7081hqe_photocathode.detect[
glass = Material('glass')
glass.set('refractive_index', 1.49)
glass.absorption_length = \
- np.array([(200, 0.1e-6), (300, 0.1e-6), (330, 1.0), (500, 2.0), (600, 1.0), (770, 0.5), (800, 0.1e-6)])
+ np.array([(200, 0.1e-6), (300, 0.1e-6), (330, 1000.0), (500, 2000.0), (600, 1000.0), (770, 500.0), (800, 0.1e-6)])
glass.set('scattering_length', 1e6)
# From WCSim
@@ -99,7 +99,7 @@ water.set('absorption_length',
11751, 8795.34, 8741.23, 7102.37, 6060.68,
4498.56, 3039.56, 2232.2, 1938, 1811.58,
1610.32, 1338.7, 1095.3, 977.525, 965.258,
- 1082.86, 876.434, 633.723, 389.87, 142.011])[::-1] / 100.0 # reversed, cm->m
+ 1082.86, 876.434, 633.723, 389.87, 142.011])[::-1] * 10.0 # reversed, cm->mm
)
water.set('scattering_length',
@@ -123,5 +123,5 @@ water.set('scattering_length',
3960.210, 3473.413, 3032.937,
2635.746, 2278.907, 1959.588,
1675.064, 1422.710, 1200.004,
- 1004.528, 833.9666, 686.1063])[::-1] / 100.0 * 0.625 # reversed, cm -> m, * magic tuning constant
+ 1004.528, 833.9666, 686.1063])[::-1] * 10.0 * 0.625 # reversed, cm -> mm, * magic tuning constant
)
diff --git a/chroma/demo/pmt.py b/chroma/demo/pmt.py
index 17af87a..9b0953f 100644
--- a/chroma/demo/pmt.py
+++ b/chroma/demo/pmt.py
@@ -5,7 +5,7 @@ from chroma.pmt import build_pmt, build_light_collector_from_file
from chroma.demo.optics import water, glass, vacuum, shiny_surface, r7081hqe_photocathode
def build_8inch_pmt(outer_material=water, nsteps=24):
- return build_pmt(dirname(__file__) + '/sno_pmt.txt', 0.003,
+ return build_pmt(dirname(__file__) + '/sno_pmt.txt', 3.0, # 3 mm
outer_material=outer_material,
glass=glass, vacuum=vacuum,
photocathode_surface=r7081hqe_photocathode,
diff --git a/chroma/generator/g4gen.py b/chroma/generator/g4gen.py
index 0ed5433..0132381 100644
--- a/chroma/generator/g4gen.py
+++ b/chroma/generator/g4gen.py
@@ -116,9 +116,9 @@ class G4Generator(object):
mass = G4ParticleTable.GetParticleTable().FindParticle(vertex.particle_name).GetPDGMass()
total_energy = vertex.ke*MeV + mass
self.particle_gun.SetParticleEnergy(total_energy)
- self.particle_gun.SetParticlePosition(G4ThreeVector(*vertex.pos)*m)
+ self.particle_gun.SetParticlePosition(G4ThreeVector(*vertex.pos)*mm)
self.particle_gun.SetParticleMomentumDirection(G4ThreeVector(*vertex.dir).unit())
- self.particle_gun.SetParticleTime(vertex.t0*s)
+ self.particle_gun.SetParticleTime(vertex.t0*ns)
if vertex.pol is not None:
self.particle_gun.SetParticlePolarization(G4ThreeVector(*vertex.pol).unit())
diff --git a/chroma/gpu/daq.py b/chroma/gpu/daq.py
index 82958dc..6ec1a79 100644
--- a/chroma/gpu/daq.py
+++ b/chroma/gpu/daq.py
@@ -20,7 +20,7 @@ class GPUChannels(object):
return event.Channels(t<1e8, t, q, self.flags.get())
class GPUDaq(object):
- def __init__(self, gpu_geometry, max_pmt_id, pmt_rms=1.2e-9):
+ def __init__(self, gpu_geometry, max_pmt_id, pmt_rms=1.2):
self.earliest_time_gpu = ga.empty(max_pmt_id+1, dtype=np.float32)
self.earliest_time_int_gpu = ga.empty(max_pmt_id+1, dtype=np.uint32)
self.channel_history_gpu = ga.zeros_like(self.earliest_time_int_gpu)
diff --git a/chroma/gpu/tools.py b/chroma/gpu/tools.py
index a11a561..38d67a1 100644
--- a/chroma/gpu/tools.py
+++ b/chroma/gpu/tools.py
@@ -2,6 +2,7 @@ import numpy as np
import pycuda.tools
from pycuda import characterize
import pycuda.driver as cuda
+import pycuda.compiler
from pycuda import gpuarray as ga
from chroma.cuda import srcdir
diff --git a/chroma/likelihood.py b/chroma/likelihood.py
index c69dc72..93a95c3 100644
--- a/chroma/likelihood.py
+++ b/chroma/likelihood.py
@@ -6,7 +6,7 @@ from chroma.tools import profile_if_possible
class Likelihood(object):
"Class to evaluate likelihoods for detector events."
- def __init__(self, sim, event=None, tbins=100, trange=(-0.5e-9, 999.5e-9),
+ def __init__(self, sim, event=None, tbins=100, trange=(-0.5, 999.5),
qbins=10, qrange=(-0.5, 9.5), time_only=True):
"""
Args:
@@ -56,7 +56,7 @@ class Likelihood(object):
hitcount, pdf_prob, pdf_prob_uncert = \
self.sim.eval_pdf(self.event.channels,
vertex_generator,
- 0.1e-9, self.trange,
+ 0.1, self.trange,
1, self.qrange,
nreps=nreps,
ndaq=ndaq,
diff --git a/chroma/pmt.py b/chroma/pmt.py
index 6ebf120..6f9e0df 100644
--- a/chroma/pmt.py
+++ b/chroma/pmt.py
@@ -34,8 +34,6 @@ def build_pmt_shell(filename, outer_material, nsteps=16):
# so that the mesh is closed
profile[0,0] = 0.0
profile[-1,0] = 0.0
- # convert mm -> m
- profile /= 1000.0
return Solid(rotate_extrude(profile[:,0], profile[:,1], nsteps), glass, outer_material, color=0xeeffffff)
@@ -52,8 +50,6 @@ def build_pmt(filename, glass_thickness, outer_material, glass,
# so that the mesh is closed
profile[0,0] = 0.0
profile[-1,0] = 0.0
- # convert mm -> m
- profile /= 1000.0
offset_profile = offset(profile, -glass_thickness)
@@ -80,9 +76,6 @@ def build_light_collector_from_file(filename, outer_material,
surface, nsteps=48):
profile = read_csv(filename)
- # Convert mm to m
- profile /= 1000.0
-
mesh = rotate_extrude(profile[:,0], profile[:,1], nsteps)
solid = Solid(mesh, outer_material, outer_material, surface=surface)
return solid
diff --git a/chroma/project.py b/chroma/project.py
index 51fba1e..05e4d2e 100644
--- a/chroma/project.py
+++ b/chroma/project.py
@@ -2,7 +2,7 @@ import numpy as np
from chroma.transform import normalize
def from_film(position=(0,0,0), axis1=(0,0,1), axis2=(1,0,0), size=(800,600),
- width=0.035, focal_length=0.018):
+ width=35.0, focal_length=18.0):
"""Project rays from a piece of film located at whose focal point is
located at `position`. `axis1` and `axis2` specify the vectors pointing
along the length and height of the film respectively.
diff --git a/src/G4chroma.cc b/src/G4chroma.cc
index 340cc4c..c0817a1 100644
--- a/src/G4chroma.cc
+++ b/src/G4chroma.cc
@@ -107,11 +107,11 @@ void PhotonTrackingAction::PreUserTrackingAction(const G4Track *track)
{
G4ParticleDefinition *particle = track->GetDefinition();
if (particle->GetParticleName() == "opticalphoton") {
- pos.push_back(track->GetPosition()/m);
+ pos.push_back(track->GetPosition()/mm);
dir.push_back(track->GetMomentumDirection());
pol.push_back(track->GetPolarization());
wavelength.push_back( (h_Planck * c_light / track->GetKineticEnergy()) / nanometer );
- t0.push_back(track->GetGlobalTime() / s);
+ t0.push_back(track->GetGlobalTime() / ns);
const_cast<G4Track *>(track)->SetTrackStatus(fStopAndKill);
}
}
diff --git a/test/data/ray_intersection.npy b/test/data/ray_intersection.npy
new file mode 100644
index 0000000..43f51c2
--- /dev/null
+++ b/test/data/ray_intersection.npy
Binary files differ
diff --git a/test/data/ray_intersection.npz b/test/data/ray_intersection.npz
deleted file mode 100644
index fba0e40..0000000
--- a/test/data/ray_intersection.npz
+++ /dev/null
Binary files differ
diff --git a/test/test_ray_intersection.py b/test/test_ray_intersection.py
index ce364e0..8a33ac3 100644
--- a/test/test_ray_intersection.py
+++ b/test/test_ray_intersection.py
@@ -9,7 +9,7 @@ class TestRayIntersection(unittest.TestCase):
self.context = chroma.gpu.create_cuda_context()
self.module = chroma.gpu.get_cu_module('mesh.h')
self.gpu_funcs = chroma.gpu.GPUFuncs(self.module)
- self.box = chroma.gpu.GPUGeometry(chroma.build(chroma.make.cube()))
+ self.box = chroma.gpu.GPUGeometry(chroma.build(chroma.make.cube(size=1000.0)))
pos, dir = chroma.project.from_film()
self.pos_gpu = ga.to_gpu(chroma.gpu.to_float3(pos))
@@ -17,10 +17,12 @@ class TestRayIntersection(unittest.TestCase):
testdir = os.path.dirname(os.path.abspath(__file__))
self.dx_standard = np.load(os.path.join(testdir,
- 'data/ray_intersection.npz'))
+ 'data/ray_intersection.npy'))
+
def test_intersection_distance(self):
dx = ga.zeros(self.pos_gpu.size, dtype=np.float32)
self.gpu_funcs.distance_to_mesh(np.int32(self.pos_gpu.size), self.pos_gpu, self.dir_gpu, self.box.gpudata, dx, block=(64,1,1), grid=(self.pos_gpu.size//64+1,1))
+
self.assertTrue((dx.get() == self.dx_standard).all())
def tearDown(self):
generated by cgit (git 2.34.1) at 2025-05-07 22:32:53 +0000