From a5193f9b6b3267ef1a3b436cdffe2dc0d3615adc Mon Sep 17 00:00:00 2001
From: Stan Seibert <stan@mtrr.org>
Date: Wed, 27 Jul 2011 16:23:27 -0400
Subject: Use the more realistic SNO database and add PMT concentrators.

---
 detectors/sno.py | 60 ++++++++++++++++++++++++++++++++++++++++++++------------
 1 file changed, 48 insertions(+), 12 deletions(-)

(limited to 'detectors/sno.py')

diff --git a/detectors/sno.py b/detectors/sno.py
index 5439739..6d33c39 100644
--- a/detectors/sno.py
+++ b/detectors/sno.py
@@ -4,9 +4,9 @@ import math
 from make import rotate_extrude
 from geometry import *
 from optics import *
-from solids import build_8inch_pmt
+from solids import build_8inch_pmt_with_lc
 from transform import rotate, make_rotation_matrix
-import sno_pmt_locations
+import os
 
 def sno_vessel(sphere_radius, neck_radius, neck_top, angle_step=math.pi/20):
     '''Compute the 2D coordinates of the profile of one side of a 
@@ -66,22 +66,58 @@ av_inside_mesh = rotate_extrude(av_inside_profile[0], av_inside_profile[1],
                                  angle_step)
 av_inside_mesh.vertices = rotate(av_inside_mesh.vertices, np.pi/2, (-1,0,0))
 
+dir = os.path.split(os.path.realpath(__file__))[0]
+
 def build_sno():
-    pmt = build_8inch_pmt()
+    import h5py
+    pmtinfo = h5py.File(dir+'/sno_pmt_info.hdf5')
+
+    pmt = build_8inch_pmt_with_lc()
+
     geo = Geometry()
-    geo.add_solid(Solid(av_outside_mesh, glass, water, color=0x99FFFFFF))
-    geo.add_solid(Solid(av_inside_mesh, water, glass, color=0x990000FF))
+    geo.add_solid(Solid(av_outside_mesh, glass, water, color=0x33FFFFFF))
+    geo.add_solid(Solid(av_inside_mesh, water, glass, color=0x330000FF))
     geo.pmtids = []
-    for x,y,z in zip(sno_pmt_locations.x, sno_pmt_locations.y,
-                     sno_pmt_locations.z):
-        direction = numpy.array((-x,-y,-z))
-        direction /= numpy.linalg.norm(direction)
+
+    snoman_id_dict = {}
+
+    pmt_offset = 122.0 - 2 * 0.01# max bucket height - 2 * facegap (mm)
+
+    for position, direction, idnum, pmt_type in zip(pmtinfo['pos'], pmtinfo['dir'],
+                                          pmtinfo['snoman_id'],
+                                          pmtinfo['pmt_type']):
+
+        # PMT type codes:
+        #        1 = normal; 2 = OWL; 3 = Low Gain; 4 = BUTT; 5 = Neck
+        #        6 = Calib channel;  10 = spare;  11 = invalid
+        if pmt_type != 1:  # All the other PMTs have nonsense directions
+            continue
+
+        
+        # Flip and renormalize
+        direction *= -1.0/numpy.linalg.norm(direction)
+
+        # Orient PMT that starts facing Y axis
         y_axis = np.array((0.0,1.0,0.0))
         axis = np.cross(direction, y_axis)
         angle = np.arccos(np.dot(y_axis, direction))
-
         rotation = make_rotation_matrix(angle, axis)
-        displacement = (x/1000.0, y/1000.0, z/1000.0)
-        geo.pmtids.append(geo.add_solid(pmt, rotation, displacement))
+
+        # Values in database are for front face of concentrator.
+        # need to shift so position is for equator of PMT
+        if pmt_type == 1:
+            displacement = position - direction * pmt_offset
+        displacement /= 1000.0 # mm -> m
+        chroma_id = geo.add_solid(pmt, rotation, displacement)
+        snoman_id_dict[chroma_id] = idnum
+        geo.pmtids.append(chroma_id)
+
+    # Convert dict to numpy array for fast array remapping later
+    chroma_ids = np.array(snoman_id_dict.keys())
+    snoman_ids = np.array(snoman_id_dict.values())
+    snoman_id_map = np.zeros(max(chroma_ids)+1) - 1 # Fill with -1 everywhere
+    snoman_id_map[chroma_ids] = snoman_ids
+
+    geo.snoman_id_map = snoman_id_map
         
     return geo
-- 
cgit