summaryrefslogtreecommitdiff
path: root/detectors/lbne.py
blob: 3c28c5849dbc38d50c5845b38b57eb140515ae49 (plain)
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 os
import sys
import numpy as np

dir = os.path.split(os.path.realpath(__file__))[0]
sys.path.append(dir + '/..')

import models
from mesh import mesh_from_stl
from solid import Solid
from geometry import Geometry
from transform import rotate, make_rotation_matrix
from itertools import product

endcap_spacing = .485

radius = 25.0/10.0
height = 50.0/10.0

nstrings = 324//10
pmts_per_string = 102//10

class LBNE(Geometry):
    def __init__(self):
        super(LBNE, self).__init__()

        pmt_mesh = mesh_from_stl(models.dir + '/hamamatsu_12inch.stl')
        pmt_mesh.vertices /= 1000.0

        pmtid = 0

        # construct the barrel
        for i in range(pmts_per_string):
            for j in range(nstrings):
                rotation = make_rotation_matrix(j*2*np.pi/nstrings, (0,0,1))
                displacement = rotate((-radius,0,-height/2+i*height/(pmts_per_string-1)), j*2*np.pi/nstrings, (0,0,1))
                self.add_solid(Solid(pmtid, pmt_mesh, rotation, displacement))
                pmtid += 1

        # construct the top endcap
        for x, y in np.array(tuple(product(\
                    np.arange(-radius+0.075, radius, endcap_spacing),
                    np.arange(-radius+0.075, radius, endcap_spacing)))):
            if np.sqrt(x**2 + y**2) <= radius:
                rotation = make_rotation_matrix(-np.pi/2, (0,1,0))
                displacement = (x,y,+height/2+height/(pmts_per_string-1)/2)
                self.add_solid(Solid(pmtid, pmt_mesh, rotation, displacement))
                pmtid += 1

        # construct the bottom endcap
        for x, y in np.array(tuple(product(\
                    np.arange(-radius+0.075, radius, endcap_spacing),
                    np.arange(-radius+0.075, radius, endcap_spacing)))):
            if np.sqrt(x**2 + y**2) <= radius:
                rotation = make_rotation_matrix(+np.pi/2, (0,1,0))
                displacement = (x,y,-height/2-height/(pmts_per_string-1)/2)
                self.add_solid(Solid(pmtid, pmt_mesh, rotation, displacement))
                pmtid += 1