From 6775d95297dae2325aaae8a706dde849d253af14 Mon Sep 17 00:00:00 2001 From: Yan-Jie Schnellbach Date: Tue, 29 Mar 2022 17:32:03 +0100 Subject: [PATCH] Added options for custom boxes, cylinders and spheres. --- sntools/detectors.py | 54 ++++++++++++++++++++++++++++++++++++++++++-- sntools/genevts.py | 19 ++++++++++++++-- 2 files changed, 69 insertions(+), 4 deletions(-) diff --git a/sntools/detectors.py b/sntools/detectors.py index 88ccada..c94a4aa 100644 --- a/sntools/detectors.py +++ b/sntools/detectors.py @@ -38,14 +38,18 @@ def wbls(x): # List of supported detector configurations supported_detectors = ["HyperK", "HyperKDR", "SuperK", "WATCHMAN", "WATCHMAN-LS", "WATCHMAN-WbLS", - "THEIA25", "THEIA100"] + "THEIA25", "THEIA100", + "CUSTOM-BOX", "CUSTOM-CYLINDER", "CUSTOM-SPHERE"] +# List of supported custom detector materials (to-do: more elegant parsing of WbLS ratios) +supported_materials = ["water", "ls", "wbls0.03", "wbls0.10"] class Detector(object): """A neutrino detector.""" - def __init__(self, name): + def __init__(self, name, arg_radius, arg_x, arg_y, arg_z, arg_material): self.name = name + if name == "HyperK": # inner detector only, 2019 optimized design self.shape = "cylinder" self.height = 6580 @@ -89,14 +93,53 @@ def __init__(self, name): self.height = 5000 self.radius = 5000 / 2 self.material = wbls(0.10) # 10% LS, 90% water + elif name == "CUSTOM-BOX": + if arg_x <= 0.0 or arg_y <= 0.0 or arg_z <= 0.0: + raise ValueError(f"Ill-defined detector measurements: x={arg_x} mm, y={arg_y} mm. z={arg_z} mm") + self.shape = "box" + # take values directly from function call + # materials will be set later for streamlining + self.x = arg_x + self.y = arg_y + self.z = arg_z + elif name == "CUSTOM-CYLINDER": + if arg_radius <= 0.0 or arg_z <= 0.0: + raise ValueError(f"Ill-defined detector measurements: r={arg_radius} mm, z={arg_z} mm") + self.shape = "cylinder" + # take values directly from function call + # materials will be set later for streamlining + self.radius = arg_radius + self.height = arg_z + elif name == "CUSTOM-SPHERE": + if arg_radius <= 0.0: + raise ValueError(f"Ill-defined detector measurements: r={arg_radius} mm") + self.shape = "sphere" + # take values directly from function call + # materials will be set later for streamlining + self.radius = arg_radius else: raise ValueError(f"Unknown detector name: {name}") + # parse material only for custom cylinder/box/sphere + if name == "CUSTOM-BOX" or name == "CUSTOM-CYLINDER" or name == "CUSTOM-SPHERE": + if arg_material == "water": + self.material = water + elif arg_material == "ls": + self.material = ls + elif arg_material == "wbls0.03": + self.material = wbls(0.03) + elif arg_material == "wbls0.10": + self.material = wbls(0.10) + else: + raise ValueError(f"Unknown custom material chosen: {arg_material}") + # calculate number of target molecules in detector if self.shape == "box": volume = self.x * self.y * self.z elif self.shape == "cylinder": volume = pi * self.radius ** 2 * self.height # assumes cylindrical detector + elif self.shape == "sphere": + volume = 4/3 * pi * self.radius ** 3 # assumes perfect sphere else: raise ValueError(f"Unknown detector shape: {self.shape}") number_density = self.material["density"] * 6.022e23 / self.material["molecular_weight"] @@ -122,6 +165,13 @@ def generate_random_vertex(self): if x ** 2 + y ** 2 < self.radius ** 2: break z = random.uniform(-self.height / 2, self.height / 2) + elif self.shape == "sphere": + while True: + x = random.uniform(-self.radius, self.radius) + y = random.uniform(-self.radius, self.radius) + z = random.uniform(-self.radius, self.radius) + if x ** 2 + y ** 2 + z ** 2 < self.radius ** 2: + break else: raise ValueError(f"Unknown detector shape: {self.shape}") return x, y, z diff --git a/sntools/genevts.py b/sntools/genevts.py index 26db277..e085492 100644 --- a/sntools/genevts.py +++ b/sntools/genevts.py @@ -18,7 +18,7 @@ import sntools from sntools.channel import gen_evts -from sntools.detectors import Detector, supported_detectors +from sntools.detectors import Detector, supported_detectors, supported_materials from sntools.formats import CompositeFlux, SNEWPYCompositeFlux from sntools.transformation import Transformation, SNEWPYTransformation @@ -99,6 +99,21 @@ def parse_command_line_options(): parser.add_argument("-d", "--detector", metavar="DETECTOR", choices=supported_detectors, default="HyperK", help="Detector configuration. Choices: %(choices)s. Default: %(default)s.") + parser.add_argument("-r", "--radius", metavar="R", type=float, default="-1", + help="Radius of custom cylindrical/spherical detector in millimeters.") + + parser.add_argument("-x", "--x", metavar="X", type=float, default="-1", + help="Width in x of custom box detector in millimeters.") + + parser.add_argument("-y", "--y", metavar="Y", type=float, default="-1", + help="Width in y of custom box detector in millimeters.") + + parser.add_argument("-z", "--z", metavar="Z", type=float, default="-1", + help="Height of custom detector in millimeters.") + + parser.add_argument("--material", metavar="MATERIAL", choices=supported_materials, default="water", + help="Custom detector material. Only used for custom detectors. Choices: %(choices)s. Default: %(default)s.") + parser.add_argument("--distance", type=float, default=10.0, help="Distance to supernova in kpc. Default: %(default)s.") parser.add_argument("--starttime", metavar="T", type=float, @@ -118,7 +133,7 @@ def parse_command_line_options(): args = parser.parse_args() - args.detector = Detector(args.detector) + args.detector = Detector(args.detector, args.radius, args.x, args.y, args.z, args.material) args.channels = args.detector.material["channel_weights"] if args.channel == "all" else [args.channel] if args.format[:7] == "SNEWPY-":