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
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
|
import numpy as np
import numpy.linalg
import math
from make import rotate_extrude
from geometry import *
from optics import *
from solids import build_8inch_pmt
from transform import rotate, make_rotation_matrix
import sno_pmt_locations
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
SNO-style acrylic vessel. The center of the sphere is at (0,0), with
the neck extending along the positive y direction.
sphere_radius: Radius of spherical part of profile
neck_radius: Radius of neck part
neck_top: y coordinate of top of neck
angle_step: angular step size (radius) between points on the sphere
Returns: Tuple of x and y coordinate numpy arrays.
'''
if neck_radius >= sphere_radius:
raise ValueError('neck_radius must be less than sphere_radius')
intersect_height = (sphere_radius**2 - neck_radius**2)**0.5
max_angle = math.atan2(intersect_height, neck_radius)
if neck_top < intersect_height:
raise ValueError('neck_top must be greater than the y-value where the sphere and cylinder intersect')
# Start with point at bottom
angles = np.arange(-math.pi/2, max_angle, angle_step)
x = list(np.cos(angles) * sphere_radius)
y = list(np.sin(angles) * sphere_radius)
x[0] = 0.0 # Round-off error might make cos(-pi/2) not exactly zero
# Neck intersection point
x.append(neck_radius)
y.append(intersect_height)
# Top of neck
x.append(neck_radius)
y.append(neck_top)
# Top of neck on axis
x.append(0.0)
y.append(neck_top)
return x, y
##### SNO Parts
angle_step = math.pi/20
av_outside_profile = sno_vessel(sphere_radius=6.0604, neck_radius=0.79375,
neck_top=10.50, angle_step=angle_step)
# For simplicity, cap the top of the AV with acrylic
av_inside_profile = sno_vessel(sphere_radius=6.0053, neck_radius=0.72898,
neck_top=10.00, angle_step=angle_step)
av_outside_mesh = rotate_extrude(av_outside_profile[0], av_outside_profile[1],
angle_step)
av_outside_mesh.vertices = rotate(av_outside_mesh.vertices, np.pi/2, (-1,0,0))
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))
def build_sno():
pmt = build_8inch_pmt()
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.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)
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))
return geo
|
