From 54102fa758a109223a9b600a9ea163bfae35a212 Mon Sep 17 00:00:00 2001 From: liviocali Date: Fri, 22 Dec 2023 21:14:09 +0100 Subject: [PATCH] Add SiPM geometry LUT --- .../light_module_desc-2.0.0.yaml | 442 ++++++++++++++++++ src/proto_nd_flow/reco/light/wvfm_sum.py | 4 +- src/proto_nd_flow/resources/geometry.py | 245 +++++++--- yamls/proto_nd_flow/resources/Geometry.yaml | 2 +- 4 files changed, 624 insertions(+), 69 deletions(-) create mode 100644 data/proto_nd_flow/light_module_desc-2.0.0.yaml diff --git a/data/proto_nd_flow/light_module_desc-2.0.0.yaml b/data/proto_nd_flow/light_module_desc-2.0.0.yaml new file mode 100644 index 00000000..a2d6c556 --- /dev/null +++ b/data/proto_nd_flow/light_module_desc-2.0.0.yaml @@ -0,0 +1,442 @@ +format_version: "0.2.0" +geometry_version: "0.2.0" + +geom: + # ArcLight + 0: { min: [-15.14 , -15.51, 0], max: [+12.03, +15.51, 0] } + # LCM + 1: { min: [-15.14 , -5.17, 0], max: [+12.03, +5.17, 0] } + +tpc_center_offset: + 0: [+15.30, 0, 0] + 1: [-15.30, 0, 0] + 2: [+15.30, 0, 0] + 3: [-15.30, 0, 0] + 4: [+15.30, 0, 0] + 5: [-15.30, 0, 0] + 6: [+15.30, 0, 0] + 7: [-15.30, 0, 0] + +det_center: + 0: [0, -46.53, -31.49] + 1: [0, -25.85, -31.49] + 2: [0, -15.51, -31.49] + 3: [0, -5.16, -31.49] + 4: [0, 15.51, -31.49] + 5: [0, 36.19, -31.49] + 6: [0, 46.53, -31.49] + 7: [0, 56.87, -31.49] + 8: [0, -46.53, +31.49] + 9: [0, -25.84, +31.49] + 10: [0, -15.51, +31.49] + 11: [0, -5.17, +31.49] + 12: [0, 15.51, +31.49] + 13: [0, 36.19, +31.49] + 14: [0, 46.53, +31.49] + 15: [0, 56.87, +31.49] + +sipm_center: + 0: [15.14, -60.07, -31.49] + 1: [15.14, -55.37, -31.49] + 2: [15.14, -48.87, -31.49] + 3: [15.14, -44.17, -31.49] + 4: [15.14, -37.67, -31.49] + 5: [15.14, -32.97, -31.49] + 6: [15.14, -29.07, -31.49] + 7: [15.14, -24.37, -31.49] + 8: [15.14, -17.87, -31.49] + 9: [15.14, -13.17, -31.49] + 10: [15.14, -6.67, -31.49] + 11: [15.14, -1.97, -31.49] + 12: [15.14, 1.97, -31.49] + 13: [15.14, 6.67, -31.49] + 14: [15.14, 13.17, -31.49] + 15: [15.14, 17.87, -31.49] + 16: [15.14, 24.37, -31.49] + 17: [15.14, 29.07, -31.49] + 18: [15.14, 32.97, -31.49] + 19: [15.14, 37.67, -31.49] + 20: [15.14, 44.17, -31.49] + 21: [15.14, 48.87, -31.49] + 22: [15.14, 55.37, -31.49] + 23: [15.14, 60.07, -31.49] + 24: [15.14, -60.07, 31.49] + 25: [15.14, -55.37, 31.49] + 26: [15.14, -48.87, 31.49] + 27: [15.14, -44.17, 31.49] + 28: [15.14, -37.67, 31.49] + 29: [15.14, -32.97, 31.49] + 30: [15.14, -29.07, 31.49] + 31: [15.14, -24.37, 31.49] + 32: [15.14, -17.87, 31.49] + 33: [15.14, -13.17, 31.49] + 34: [15.14, -6.67, 31.49] + 35: [15.14, -1.97, 31.49] + 36: [15.14, 1.97, 31.49] + 37: [15.14, 6.67, 31.49] + 38: [15.14, 13.17, 31.49] + 39: [15.14, 17.87, 31.49] + 40: [15.14, 24.37, 31.49] + 41: [15.14, 29.07, 31.49] + 42: [15.14, 32.97, 31.49] + 43: [15.14, 37.67, 31.49] + 44: [15.14, 44.17, 31.49] + 45: [15.14, 48.87, 31.49] + 46: [15.14, 55.37, 31.49] + 47: [15.14, 60.07, 31.49] + +ch_to_vert_bin: # 0:ACL 1:LCM + 0: [-1,-1,-1,-1,0,1,2,3,4,5,12,13,14,15,16,17,-1,-1,-1,-1,0,1,2,3,4,5,12,13,14,15,16,17, -1,-1,-1,-1,0,1,2,3,4,5,12,13,14,15,16,17,-1,-1,-1,-1,0,1,2,3,4,5,12,13,14,15,16,17] + 1: [-1,-1,-1,-1,6,7,8,9,10,11,18,19,20,21,22,23,-1,-1,-1,-1,6,7,8,9,10,11,18,19,20,21,22,23,-1,-1,-1,-1,6,7,8,9,10,11,18,19,20,21,22,23,-1,-1,-1,-1,6,7,8,9,10,11,18,19,20,21,22,23] + +adc_to_det_type: + 0: 0 + 1: 1 + 2: 0 + 3: 1 + 4: 0 + 5: 1 + 6: 0 + 7: 1 + +det_side: #TPC side 0: -z direction 1: +z direction + 0: 0 + 1: 0 + 2: 0 + 3: 0 + 4: 0 + 5: 0 + 6: 0 + 7: 0 + 8: 1 + 9: 1 + 10: 1 + 11: 1 + 12: 1 + 13: 1 + 14: 1 + 15: 1 + +det_geom: + 0: 0 + 1: 1 + 2: 1 + 3: 1 + 4: 0 + 5: 1 + 6: 1 + 7: 1 + 8: 0 + 9: 1 + 10: 1 + 11: 1 + 12: 0 + 13: 1 + 14: 1 + 15: 1 + +det_adc: + 0: #TPC + # -z, increasing y (det-> adc) + 0: 0 + 1: 1 + 2: 1 + 3: 1 + 4: 0 + 5: 1 + 6: 1 + 7: 1 + # +z, increasing y + 8: 0 + 9: 1 + 10: 1 + 11: 1 + 12: 0 + 13: 1 + 14: 1 + 15: 1 + + 1: + # -z, increasing y + 0: 0 + 1: 1 + 2: 1 + 3: 1 + 4: 0 + 5: 1 + 6: 1 + 7: 1 + # +z, increasing y + 8: 0 + 9: 1 + 10: 1 + 11: 1 + 12: 0 + 13: 1 + 14: 1 + 15: 1 + + 2: + # -z, increasing y + 0: 2 + 1: 3 + 2: 3 + 3: 3 + 4: 2 + 5: 3 + 6: 3 + 7: 3 + # +z, increasing y + 8: 2 + 9: 3 + 10: 3 + 11: 3 + 12: 2 + 13: 3 + 14: 3 + 15: 3 + + 3: + # -z, increasing y + 0: 2 + 1: 3 + 2: 3 + 3: 3 + 4: 2 + 5: 3 + 6: 3 + 7: 3 + # +z, increasing y + 8: 2 + 9: 3 + 10: 3 + 11: 3 + 12: 2 + 13: 3 + 14: 3 + 15: 3 + + 4: + # -z, increasing y + 0: 4 + 1: 5 + 2: 5 + 3: 5 + 4: 4 + 5: 5 + 6: 5 + 7: 5 + # +z, increasing y + 8: 4 + 9: 5 + 10: 5 + 11: 5 + 12: 4 + 13: 5 + 14: 5 + 15: 5 + + 5: + # -z, increasing y + 0: 4 + 1: 5 + 2: 5 + 3: 5 + 4: 4 + 5: 5 + 6: 5 + 7: 5 + # +z, increasing y + 8: 4 + 9: 5 + 10: 5 + 11: 5 + 12: 4 + 13: 5 + 14: 5 + 15: 5 + + 6: + # -z, increasing y + 0: 6 + 1: 7 + 2: 7 + 3: 7 + 4: 6 + 5: 7 + 6: 7 + 7: 7 + # +z, increasing y + 8: 6 + 9: 7 + 10: 7 + 11: 7 + 12: 6 + 13: 7 + 14: 7 + 15: 7 + + 7: + # -z, increasing y + 0: 6 + 1: 7 + 2: 7 + 3: 7 + 4: 6 + 5: 7 + 6: 7 + 7: 7 + # +z, increasing y + 8: 6 + 9: 7 + 10: 7 + 11: 7 + 12: 6 + 13: 7 + 14: 7 + 15: 7 + +det_chan: + 0: + 0: [4,5,6,7,8,9] + 1: [4,5] + 2: [6,7] + 3: [8,9] + 4: [10,11,12,13,14,15] + 5: [10,11] + 6: [12,13] + 7: [14,15] + 8: [20,21,22,23,24,25] + 9: [20,21] + 10: [22,23] + 11: [24,25] + 12: [26,27,28,29,30,31] + 13: [26,27] + 14: [28,29] + 15: [30,31] + + 1: + 0: [52,53,54,55,56,57] + 1: [52,53] + 2: [54,55] + 3: [56,57] + 4: [58,59,60,61,62,63] + 5: [58,59] + 6: [60,61] + 7: [62,63] + 8: [36,37,38,39,40,41] + 9: [36,37] + 10: [38,39] + 11: [40,41] + 12: [42,43,44,45,46,47] + 13: [42,43] + 14: [44,45] + 15: [46,47] + + 2: + 0: [4,5,6,7,8,9] + 1: [4,5] + 2: [6,7] + 3: [8,9] + 4: [10,11,12,13,14,15] + 5: [10,11] + 6: [12,13] + 7: [14,15] + 8: [20,21,22,23,24,25] + 9: [20,21] + 10: [22,23] + 11: [24,25] + 12: [26,27,28,29,30,31] + 13: [26,27] + 14: [28,29] + 15: [30,31] + + 3: + 0: [52,53,54,55,56,57] + 1: [52,53] + 2: [54,55] + 3: [56,57] + 4: [58,59,60,61,62,63] + 5: [58,59] + 6: [60,61] + 7: [62,63] + 8: [36,37,38,39,40,41] + 9: [36,37] + 10: [38,39] + 11: [40,41] + 12: [42,43,44,45,46,47] + 13: [42,43] + 14: [44,45] + 15: [46,47] + + 4: + 0: [4,5,6,7,8,9] + 1: [4,5] + 2: [6,7] + 3: [8,9] + 4: [10,11,12,13,14,15] + 5: [10,11] + 6: [12,13] + 7: [14,15] + 8: [20,21,22,23,24,25] + 9: [20,21] + 10: [22,23] + 11: [24,25] + 12: [26,27,28,29,30,31] + 13: [26,27] + 14: [28,29] + 15: [30,31] + + 5: + 0: [52,53,54,55,56,57] + 1: [52,53] + 2: [54,55] + 3: [56,57] + 4: [58,59,60,61,62,63] + 5: [58,59] + 6: [60,61] + 7: [62,63] + 8: [36,37,38,39,40,41] + 9: [36,37] + 10: [38,39] + 11: [40,41] + 12: [42,43,44,45,46,47] + 13: [42,43] + 14: [44,45] + 15: [46,47] + + 6: + 0: [4,5,6,7,8,9] + 1: [4,5] + 2: [6,7] + 3: [8,9] + 4: [10,11,12,13,14,15] + 5: [10,11] + 6: [12,13] + 7: [14,15] + 8: [20,21,22,23,24,25] + 9: [20,21] + 10: [22,23] + 11: [24,25] + 12: [26,27,28,29,30,31] + 13: [26,27] + 14: [28,29] + 15: [30,31] + + 7: + 0: [52,53,54,55,56,57] + 1: [52,53] + 2: [54,55] + 3: [56,57] + 4: [58,59,60,61,62,63] + 5: [58,59] + 6: [60,61] + 7: [62,63] + 8: [36,37,38,39,40,41] + 9: [36,37] + 10: [38,39] + 11: [40,41] + 12: [42,43,44,45,46,47] + 13: [42,43] + 14: [44,45] + 15: [46,47] \ No newline at end of file diff --git a/src/proto_nd_flow/reco/light/wvfm_sum.py b/src/proto_nd_flow/reco/light/wvfm_sum.py index e21e3fbf..4d4afaf5 100644 --- a/src/proto_nd_flow/reco/light/wvfm_sum.py +++ b/src/proto_nd_flow/reco/light/wvfm_sum.py @@ -87,7 +87,7 @@ def run(self, source_name, source_slice, cache): for adc in range(wvfm_data['samples'].shape[1]): for chan in range(wvfm_data['samples'].shape[2]): - tpc_id = resources['Geometry'].tpc_id[(adc,chan)] + tpc_id = resources['Geometry'].sipm_rel_pos[(adc,chan)][0][0] det_id = resources['Geometry'].det_id[(adc,chan)] if tpc_id < 0 or det_id < 0: continue @@ -98,7 +98,7 @@ def run(self, source_name, source_slice, cache): for adc in range(wvfm_data['samples'].shape[1]): for chan in range(wvfm_data['samples'].shape[2]): - tpc_id = resources['Geometry'].tpc_id[(adc,chan)] + tpc_id = resources['Geometry'].sipm_rel_pos[(adc,chan)][0][0] det_id = resources['Geometry'].det_id[(adc,chan)] if tpc_id < 0 or det_id < 0: continue diff --git a/src/proto_nd_flow/resources/geometry.py b/src/proto_nd_flow/resources/geometry.py index 046b9ee5..53df3cbe 100644 --- a/src/proto_nd_flow/resources/geometry.py +++ b/src/proto_nd_flow/resources/geometry.py @@ -53,9 +53,11 @@ class Geometry(H5FlowResource): is in the LAr fiducial volume Provides (for light geometry): - - ``tpc_id``: lookup table for TPC number for light detectors + - ``det_rel_pos``: lookup table for relative position (TPC,side,vertical position from bottom) of light detectors (Full ArCLight or LCM) + - ``sipm_rel_pos``: lookup table for lookup table for relative position (TPC,side,vertical position from bottom) of SiPMs (Single SiPM) - ``det_id``: lookup table for detector number from adc, channel id - ``det_bounds``: lookup table for detector minimum and maximum corners light detectors + - ``sipm_abs_pos``: lookup table for sipm absolute position (x,y,z) in cm - ``solid_angle()``: helper function for determining the solid angle of a given detector Example usage:: @@ -113,6 +115,15 @@ def init(self, source_name): if not self.data: # first time loading geometry, save to file + with open(self.crs_geometry_file) as gf: + self.crs_geometry_yaml = yaml.load(gf, Loader=yaml.FullLoader) + + with open(self.det_geometry_file) as dgf: + self.det_geometry_yaml = yaml.load(dgf, Loader=yaml.FullLoader) + + with open(self.lrs_geometry_file) as gf: + self.lrs_geometry_yaml = yaml.load(gf, Loader=yaml.FullLoader) + self.load_geometry() self.data_manager.set_attrs(self.path, @@ -134,9 +145,11 @@ def init(self, source_name): write_lut(self.data_manager, self.path, self.pixel_coordinates_2D, 'pixel_coordinates_2D') write_lut(self.data_manager, self.path, self.tile_id, 'tile_id') - write_lut(self.data_manager, self.path, self.tpc_id, 'tpc_id') + write_lut(self.data_manager, self.path, self.det_rel_pos, 'det_rel_pos') + write_lut(self.data_manager, self.path, self.sipm_rel_pos, 'sipm_rel_pos') write_lut(self.data_manager, self.path, self.det_id, 'det_id') write_lut(self.data_manager, self.path, self.det_bounds, 'det_bounds') + write_lut(self.data_manager, self.path, self.sipm_abs_pos, 'sipm_abs_pos') else: assert_compat_version(self.class_version, self.data['class_version']) @@ -152,14 +165,17 @@ def init(self, source_name): self._pixel_coordinates_2D = read_lut(self.data_manager, self.path, 'pixel_coordinates_2D') self._tile_id = read_lut(self.data_manager, self.path, 'tile_id') - self._tpc_id = read_lut(self.data_manager, self.path, 'tpc_id') + self._det_rel_pos = read_lut(self.data_manager, self.path, 'det_rel_pos') + self._sipm_rel_pos = read_lut(self.data_manager, self.path, 'sipm_rel_pos') self._det_id = read_lut(self.data_manager, self.path, 'det_id') self._det_bounds = read_lut(self.data_manager, self.path, 'det_bounds') + self._sipm_abs_pos = read_lut(self.data_manager, self.path, 'sipm_abs_pos') lut_size = (self.anode_drift_coordinate.nbytes + self.drift_dir.nbytes + self.pixel_coordinates_2D.nbytes + self.tile_id.nbytes - + self.tpc_id.nbytes + self.det_id.nbytes - + self.det_bounds.nbytes) + + self.det_rel_pos.nbytes + self.det_rel_pos.nbytes + + self.det_id.nbytes + self.det_bounds.nbytes + + self.sipm_abs_pos.nbytes) if self.rank == 0: logging.info(f'Geometry LUT(s) size: {lut_size/1024/1024:0.02f}MB') @@ -328,10 +344,8 @@ def _get_module_RO_bounds(self): ''' Get module_RO_bounds from pre-saved 2D pixel coordinates and anode drift coordinates ''' - with open(self.det_geometry_file) as dgf: - det_geometry_yaml = yaml.load(dgf, Loader=yaml.FullLoader) - module_to_io_groups = det_geometry_yaml['module_to_io_groups'] + module_to_io_groups = self.det_geometry_yaml['module_to_io_groups'] self._module_RO_bounds = [] @@ -393,20 +407,30 @@ def _rotate_pixel(pixel_pos, tile_orientation): ## Light geometry methods ## @property - def tpc_id(self): + def det_rel_pos(self): ''' - Lookup table for TPC id, usage:: + Lookup table for detector relative position, usage:: - resource['Geometry'].tpc_id[(adc_index, channel_index)] + resource['Geometry'].det_rel_pos[(tpc_index, detector_index)] ''' - return self._tpc_id + return self._det_rel_pos + + @property + def sipm_rel_pos(self): + ''' + Lookup table for detector relative position, usage:: + + resource['Geometry'].sipm_rel_pos[(adc_index, channel_index)] + + ''' + return self._det_rel_pos @property def det_id(self): ''' - Lookup table for detector id within a TPC, usage:: + Lookup table for TPC and detector id, usage:: resource['Geometry'].det_id[(adc_index, channel_index)] @@ -423,6 +447,16 @@ def det_bounds(self): ''' return self._det_bounds + + @property + def sipm_abs_pos(self): + ''' + Lookup table for SiPM center xyz coordinate, usage:: + + resource['Geometry'].sipm_abs_pos[(adc_index, channel_index)] + + ''' + return self._sipm_abs_pos @staticmethod @@ -472,6 +506,72 @@ def solid_angle(self, xyz, tpc_id, det_id): omega += det_y_sign * det_z_sign * np.arctan2(np.abs(det_y-y) * np.abs(det_z-z), np.abs(det_x-x)* d) return omega + def get_sipm_rel_pos(self, adc, channel): + ''' Returns + - TPC number starting at 0 (Attention not to be confused wiht CRS IO group) + - TPC side with 0: -z direction 1: +z direction + - Vertical position starting from bottom + + if channel not used, returns NaN,NaN,NaN + ''' + + # Get TPC/det number + tpc = -1 + det = -1 + for tpc_temp, det_map in self.lrs_geometry_yaml["det_adc"].items(): + for det_temp, adc_map in det_map.items(): + if adc_map == adc: + if channel in self.lrs_geometry_yaml["det_chan"][tpc_temp][det_temp]: + tpc = tpc_temp + det = det_temp + + # Return NaN if adc-channel does not exist + if tpc == -1 or det == -1: + return [-1,-1,-1] + + det_type = self.lrs_geometry_yaml["adc_to_det_type"][adc] + + # Get TPC side + side = self.lrs_geometry_yaml["det_side"][det] + + # Get vertical position + # Get Y pos + if det_type == 0: + vert_pos = self.lrs_geometry_yaml["ch_to_vert_bin"][0][channel] + else: + vert_pos = self.lrs_geometry_yaml["ch_to_vert_bin"][1][channel] + + return tpc, side, vert_pos + + def get_sipm_abs_pos(self,adc,channel): + '''Returns x,y,z position of each SiPM in cm (z=beam direction) + if channel not used, returns NaN,NaN,NaN''' + + tpc, side, vert_pos = self.get_sipm_rel_pos(adc,channel) + tpc_channel = vert_pos + side*(len(self.lrs_geometry_yaml["sipm_center"])//2) + + if np.isnan(tpc): + return [-1,-1,-1] + + # Get X pos + x_pos = self.det_geometry_yaml["tpc_offsets"][tpc//2][0] + self.lrs_geometry_yaml["tpc_center_offset"][tpc][0] + if tpc % 2 == 0: + x_pos += self.lrs_geometry_yaml["sipm_center"][tpc_channel][0] + else: + x_pos -= self.lrs_geometry_yaml["sipm_center"][tpc_channel][0] + + # Get Y pos + y_pos = self.det_geometry_yaml["tpc_offsets"][tpc//2][1] + self.lrs_geometry_yaml["tpc_center_offset"][tpc][1] + y_pos += self.lrs_geometry_yaml["sipm_center"][tpc_channel][1] + + # Get Z pos + z_pos = self.det_geometry_yaml["tpc_offsets"][tpc//2][2] + self.lrs_geometry_yaml["tpc_center_offset"][tpc][2] + if tpc % 2 == 0: + z_pos += self.lrs_geometry_yaml["sipm_center"][tpc_channel][2] + else: + z_pos -= self.lrs_geometry_yaml["sipm_center"][tpc_channel][2] + + return x_pos, y_pos, z_pos ## Load light and charge geometry ## @@ -484,79 +584,92 @@ def _load_light_geometry(self): if self.rank == 0: logging.warning(f'Loading geometry from {self.lrs_geometry_file}...') - with open(self.lrs_geometry_file) as gf: - geometry = yaml.load(gf, Loader=yaml.FullLoader) + assert_compat_version(self.lrs_geometry_yaml['format_version'], '0.2.0') + + mod_ids = np.array([v for v in self.det_geometry_yaml['module_to_tpcs'].keys()]) + tpc_ids = np.array([v for v in self.lrs_geometry_yaml['tpc_center_offset'].keys()]) + det_ids = np.array([v for v in self.lrs_geometry_yaml['det_center'].keys()]) + adc_ids = np.array([v for v in self.lrs_geometry_yaml['adc_to_det_type'].keys()]) + max_chan_per_det = max([len(chan) for tpc in self.lrs_geometry_yaml['det_chan'].values() for chan in tpc.values()]) + chan_ids = np.unique(sum([chan for tpc in self.lrs_geometry_yaml['det_chan'].values() for chan in tpc.values()],[])) - # enforce that light geometry formatting is as expected - assert_compat_version(geometry['format_version'], '0.0.0') + tpc_mod = np.full(tpc_ids.shape, -1, dtype=int) + for i, mod in enumerate(self.det_geometry_yaml["module_to_tpcs"]): + for j, tpc in enumerate(self.det_geometry_yaml["module_to_tpcs"][mod]): + tpc_mod[tpc] = i - tpc_ids = np.array([v for v in geometry['tpc_center'].keys()]) - det_ids = np.array([v for v in geometry['det_center'].keys()]) - max_chan = max([len(chan) for tpc in geometry['det_chan'].values() for chan in tpc.values()]) + det_min_max = [(min(tpc_ids), max(tpc_ids)), + (min(det_ids), max(det_ids))] + self._det_rel_pos = LUT('i4', *det_min_max, shape=(3,)) + self._det_rel_pos.default = -1 shape = tpc_ids.shape + det_ids.shape det_adc = np.full(shape, -1, dtype=int) - det_chan = np.full(shape + (max_chan,), -1, dtype=int) - det_chan_mask = np.zeros(shape + (max_chan,), dtype=bool) + det_side = np.full(shape, -1, dtype=int) + det_vert_pos = np.full(shape, -1, dtype=int) + det_chan = np.full(shape + (max_chan_per_det,), -1, dtype=int) + det_chan_mask = np.zeros(shape + (max_chan_per_det,), dtype=bool) det_bounds = np.zeros(shape + (2,3), dtype=float) for i, tpc in enumerate(tpc_ids): for j, det in enumerate(det_ids): - det_adc[i,j] = geometry['det_adc'][tpc][det] - det_chan[i,j,:len(geometry['det_chan'][tpc][det])] = geometry['det_chan'][tpc][det] - - tpc_center = np.array(geometry['tpc_center'][tpc]) - det_geom = geometry['geom'][geometry['det_geom'][det]] - det_center = np.array(geometry['det_center'][det]) + det_adc[i,j] = self.lrs_geometry_yaml['det_adc'][tpc][det] + det_side[i,j] = self.lrs_geometry_yaml['det_side'][det] + det_vert_pos[i,j] = [key for key, value in self.lrs_geometry_yaml['det_side'].items() if value == det_side[i,j]].index(det) + det_chan[i,j,:len(self.lrs_geometry_yaml['det_chan'][tpc][det])] = self.lrs_geometry_yaml['det_chan'][tpc][det] + + tpc_center = (np.array(self.lrs_geometry_yaml['tpc_center_offset'][tpc]) + + np.array(self.det_geometry_yaml["tpc_offsets"][tpc_mod[i]])) + det_geom = self.lrs_geometry_yaml['geom'][self.lrs_geometry_yaml['det_geom'][det]] + det_center = np.array(self.lrs_geometry_yaml['det_center'][det]) det_bounds[i,j,0] = tpc_center + det_center + np.array(det_geom['min']) det_bounds[i,j,1] = tpc_center + det_center + np.array(det_geom['max']) + self._det_rel_pos[i,j] = np.array((tpc,det_side[i,j],det_vert_pos[i,j])) det_chan_mask = det_chan != -1 det_adc, det_chan, tpc_ids, det_ids = np.broadcast_arrays( - det_adc[...,np.newaxis], det_chan, - tpc_ids[...,np.newaxis,np.newaxis], det_ids[...,np.newaxis]) + det_adc[...,np.newaxis], + det_chan, tpc_ids[...,np.newaxis,np.newaxis], det_ids[...,np.newaxis]) + + adc_chan_min_max = [(min(adc_ids), max(adc_ids)), + (min(chan_ids), max(chan_ids))] + self._sipm_abs_pos = LUT('f4', *adc_chan_min_max, shape=(3,)) + self._sipm_abs_pos.default = -1 - adc_chan_min_max = [(min(det_adc[det_chan_mask]), max(det_adc[det_chan_mask])), - (min(det_chan[det_chan_mask]), max(det_chan[det_chan_mask]))] - self._tpc_id = LUT('i4', *adc_chan_min_max) - self._tpc_id.default = -1 + self._sipm_rel_pos = LUT('i4', *adc_chan_min_max, shape=(3,)) + self._sipm_rel_pos.default = -1 self._det_id = LUT('i4', *adc_chan_min_max) self._det_id.default = -1 - det_min_max = [(min(tpc_ids[det_chan_mask]), max(tpc_ids[det_chan_mask])), - (min(det_ids[det_chan_mask]), max(det_ids[det_chan_mask]))] self._det_bounds = LUT('f4', *det_min_max, shape=(2,3)) self._det_bounds.default = 0. - self._tpc_id[(det_adc[det_chan_mask], det_chan[det_chan_mask])] = tpc_ids[det_chan_mask] self._det_id[(det_adc[det_chan_mask], det_chan[det_chan_mask])] = det_ids[det_chan_mask] + for adc in adc_ids: + for chan in chan_ids: + self._sipm_rel_pos = np.array(self.get_sipm_rel_pos) + self._sipm_abs_pos[(adc,chan)] = np.array(self.get_sipm_abs_pos(adc,chan)) + tpc_ids, det_ids, det_chan_mask = tpc_ids[...,0], det_ids[...,0], det_chan_mask[...,0] self._det_bounds[(tpc_ids[det_chan_mask], det_ids[det_chan_mask])] = det_bounds[det_chan_mask] - def _load_charge_geometry(self): if self.rank == 0: logging.warning(f'Loading geometry from {self.crs_geometry_file}...') - - with open(self.crs_geometry_file) as gf: - geometry_yaml = yaml.load(gf, Loader=yaml.FullLoader) - - with open(self.det_geometry_file) as dgf: - det_geometry_yaml = yaml.load(dgf, Loader=yaml.FullLoader) - - if 'multitile_layout_version' not in geometry_yaml.keys(): + + if 'multitile_layout_version' not in self.crs_geometry_yaml.keys(): raise RuntimeError('Only multi-tile geometry configurations are accepted') - self._pixel_pitch = geometry_yaml['pixel_pitch'] / units.cm # convert mm -> cm - self._max_drift_distance = det_geometry_yaml['drift_length'] # det geo yaml is already in cm - chip_channel_to_position = geometry_yaml['chip_channel_to_position'] - tile_orientations = geometry_yaml['tile_orientations'] - tile_positions = geometry_yaml['tile_positions'] - mod_centers = det_geometry_yaml['tpc_offsets'] - tile_chip_to_io = geometry_yaml['tile_chip_to_io'] - module_to_io_groups = det_geometry_yaml['module_to_io_groups'] + self._pixel_pitch = self.crs_geometry_yaml['pixel_pitch'] / units.cm # convert mm -> cm + self._max_drift_distance = self.det_geometry_yaml['drift_length'] # det geo yaml is already in cm + chip_channel_to_position = self.crs_geometry_yaml['chip_channel_to_position'] + tile_orientations = self.crs_geometry_yaml['tile_orientations'] + tile_positions = self.crs_geometry_yaml['tile_positions'] + mod_centers = self.det_geometry_yaml['tpc_offsets'] + tile_chip_to_io = self.crs_geometry_yaml['tile_chip_to_io'] + module_to_io_groups = self.det_geometry_yaml['module_to_io_groups'] zs = np.array(list(chip_channel_to_position.values()))[:, 0] * self.pixel_pitch ys = np.array(list(chip_channel_to_position.values()))[:, 1] * self.pixel_pitch @@ -565,27 +678,27 @@ def _load_charge_geometry(self): tile_geometry = {} - tiles = np.arange(1,len(geometry_yaml['tile_chip_to_io'])*len(det_geometry_yaml['module_to_io_groups'])+1) + tiles = np.arange(1,len(self.crs_geometry_yaml['tile_chip_to_io'])*len(self.det_geometry_yaml['module_to_io_groups'])+1) io_groups = [ - geometry_yaml['tile_chip_to_io'][tile][chip] // 1000 + (mod-1)*2 - for tile in geometry_yaml['tile_chip_to_io'] - for chip in geometry_yaml['tile_chip_to_io'][tile] + self.crs_geometry_yaml['tile_chip_to_io'][tile][chip] // 1000 + (mod-1)*2 + for tile in self.crs_geometry_yaml['tile_chip_to_io'] + for chip in self.crs_geometry_yaml['tile_chip_to_io'][tile] for mod in module_to_io_groups ] io_channels = [ - geometry_yaml['tile_chip_to_io'][tile][chip] % 1000 - for tile in geometry_yaml['tile_chip_to_io'] - for chip in geometry_yaml['tile_chip_to_io'][tile] + self.crs_geometry_yaml['tile_chip_to_io'][tile][chip] % 1000 + for tile in self.crs_geometry_yaml['tile_chip_to_io'] + for chip in self.crs_geometry_yaml['tile_chip_to_io'][tile] for mod in module_to_io_groups ] chip_ids = [ chip_channel // 1000 - for chip_channel in geometry_yaml['chip_channel_to_position'] + for chip_channel in self.crs_geometry_yaml['chip_channel_to_position'] for mod in module_to_io_groups ] channel_ids = [ chip_channel % 1000 - for chip_channel in geometry_yaml['chip_channel_to_position'] + for chip_channel in self.crs_geometry_yaml['chip_channel_to_position'] for mod in module_to_io_groups ] @@ -617,7 +730,7 @@ def _load_charge_geometry(self): for chip in tile_chip_to_io[tile]: io_group_io_channel = tile_chip_to_io[tile][chip] - io_group = io_group_io_channel//1000 + (module_id-1)*len(det_geometry_yaml['module_to_io_groups'][module_id]) + io_group = io_group_io_channel//1000 + (module_id-1)*len(self.det_geometry_yaml['module_to_io_groups'][module_id]) io_channel = io_group_io_channel % 1000 self._tile_id[([io_group], [io_channel])] = tile+(module_id-1)*len(tile_chip_to_io) @@ -638,11 +751,11 @@ def _load_charge_geometry(self): if self.network_agnostic == True: warnings.warn('Encountered an out-of-network chip, but because you enabled ``network_agnostic``, we will carry on with assumptions about the io group and io channel') # using the info about the first chip on the tile for all others - io_group_io_channel = list(geometry_yaml['tile_chip_to_io'][tile].values())[0] + io_group_io_channel = list(self.crs_geometry_yaml['tile_chip_to_io'][tile].values())[0] else: continue - io_group = io_group_io_channel // 1000 + (module_id-1)*len(det_geometry_yaml['module_to_io_groups'][module_id]) + io_group = io_group_io_channel // 1000 + (module_id-1)*len(self.det_geometry_yaml['module_to_io_groups'][module_id]) io_channel = io_group_io_channel % 1000 z = chip_channel_to_position[chip_channel][0] * \ diff --git a/yamls/proto_nd_flow/resources/Geometry.yaml b/yamls/proto_nd_flow/resources/Geometry.yaml index 94c72587..f57ce831 100644 --- a/yamls/proto_nd_flow/resources/Geometry.yaml +++ b/yamls/proto_nd_flow/resources/Geometry.yaml @@ -4,6 +4,6 @@ params: path: 'geometry_info' det_geometry_file: 'data/proto_nd_flow/2x2.yaml' crs_geometry_file: 'data/proto_nd_flow/multi_tile_layout-2.4.16.yaml' - lrs_geometry_file: 'data/proto_nd_flow/light_module_desc-1.0.0.yaml' + lrs_geometry_file: 'data/proto_nd_flow/light_module_desc-2.0.0.yaml' beam_direction: 'z' drift_direction: 'x'