Added options for custom boxes, cylinders and spheres.

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

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-":