updated S2 corrected area

straxen/corrections_services.py

@@ -14,10 +14,10 @@
 
 export, __all__ = strax.exporter()
 
-corrections_w_file = ['mlp_model', 'gcn_model', 'cnn_model',
-                      's2_xy_map', 's1_xyz_map_mlp', 's1_xyz_map_cnn',
-                      's1_xyz_map_gcn', 'fdc_map_mlp', 'fdc_map_gcn',
-                      'fdc_map_cnn']
+corrections_w_file = ['mlp_model', 'cnn_model', 'gcn_model',
+                      's2_xy_map_mlp', 's2_xy_map_cnn', 's2_xy_map_gcn', 's2_xy_map',
+                      's1_xyz_map_mlp', 's1_xyz_map_cnn', 's1_xyz_map_gcn',
+                      'fdc_map_mlp', 'fdc_map_cnn', 'fdc_map_gcn']
 
 single_value_corrections = ['elife_xenon1t', 'elife', 'baseline_samples_nv',
                             'electron_drift_velocity', 'electron_drift_time_gate']
@@ -27,7 +27,7 @@
 
 # needed because we pass these names as strax options which then get paired with the default reconstruction algorithm
 # important for apply_cmt_version
-posrec_corrections_basenames = ['s1_xyz_map', 'fdc_map']
+posrec_corrections_basenames = ['s1_xyz_map', 'fdc_map', 's2_xy_map']
 
 
 class CMTVersionError(Exception):

straxen/plugins/event_processing.py

@@ -599,20 +599,20 @@ def get_veto_tags(events, split_tags, result):
             (first_sr1_run, pax_file('XENON1T_s1_xyz_lce_true_kr83m_SR1_pax-680_fdc-3d_v0.json'))]),  # noqa
     strax.Option(
         's2_xy_correction_map',
-        help="S2 (x, y) correction map. Correct S2 position dependence "
+        help="S2 (x, y) correction map. Correct S2 position dependence, including S2 top, bottom, and total."
              "manly due to bending of anode/gate-grid, PMT quantum efficiency "
              "and extraction field distribution, as well as other geometric factors.",
         default_by_run=[
             (0, pax_file('XENON1T_s2_xy_ly_SR0_24Feb2017.json')),
             (170118_1327, pax_file('XENON1T_s2_xy_ly_SR1_v2.2.json'))]),
-   strax.Option(
+    strax.Option(
         'elife_conf',
         default=("elife", "ONLINE", True),
         help='Electron lifetime '
              'Specify as (model_type->str, model_config->str, is_nT->bool) '
              'where model_type can be "elife" or "elife_constant" '
              'and model_config can be a version.'
-   ),
+    ),
     *DEFAULT_POSREC_ALGO_OPTION
 )
 class CorrectedAreas(strax.Plugin):
@@ -626,14 +626,27 @@ class CorrectedAreas(strax.Plugin):
     Note:
         Please be aware that for both, the main and alternative S1, the
         area is corrected according to the xy-position of the main S2.
+
+        There are now 3 components of cS2s: cs2_top, cS2_bottom and cs2.
+        cs2_top and cs2_bottom are corrected by the corresponding maps,
+        and cs2 is the sum of the two.
     """
-    __version__ = '0.1.1'
+    __version__ = '0.1.2'
 
     depends_on = ['event_basics', 'event_positions']
     dtype = [('cs1', np.float32, 'Corrected S1 area [PE]'),
              ('cs2', np.float32, 'Corrected S2 area [PE]'),
+             ('cs2_top', np.float32, 'Corrected S2 area in the top PMT array [PE]'),
+             ('cs2_bottom', np.float32, 'Corrected S2 area in the bottom PMT array [PE]'),
              ('alt_cs1', np.float32, 'Corrected area of the alternate S1 [PE]'),
-             ('alt_cs2', np.float32, 'Corrected area of the alternate S2 [PE]')
+             ('alt_cs2', np.float32, 'Corrected area of the alternate S2 [PE]'),
+             ('alt_cs2_top', np.float32, 'Corrected area of the alternate S2 in the top PMT array [PE]'),
+             ('alt_cs2_bottom', np.float32, 'Corrected area of the alternate S2 in the bottom PMT array [PE]'),
+             ('cs2_area_fraction_top', np.float32, 'Main cS2 fraction of area seen by the top PMT array'),
+             ('alt_cs2_area_fraction_top', np.float32, 'Alternate cS2 fraction of area seen by the top PMT array'),
+             ('elife_correction', np.float32, 'Main cS2 correction factor due to electron lifetime'),
+             ('alt_elife_correction', np.float32, 'Alt cS2 correction factor due to electron lifetime'),
+
              ] + strax.time_fields
 
     def setup(self):
@@ -653,37 +666,54 @@ def setup(self):
 
         self.s2_map = InterpolatingMap(
             get_cmt_resource(self.run_id,
-                             tuple([*self.config['s2_xy_correction_map']]),
-                             fmt='text'))
+                             tuple(['suffix',
+                                    self.config['default_reconstruction_algorithm'],
+                                    *self.config['s2_xy_correction_map']]),
+                             fmt='json'))
 
     def compute(self, events):
+        result = dict(
+            time=events['time'],
+            endtime=strax.endtime(events)
+        )
+
         # S1 corrections depend on the actual corrected event position.
         # We use this also for the alternate S1; for e.g. Kr this is
         # fine as the S1 correction varies slowly.
         event_positions = np.vstack([events['x'], events['y'], events['z']]).T
+        for peak_type in ["", "alt_"]:
+            result[f"{peak_type}cs1"] = events[f'{peak_type}s1_area'] / self.s1_map(event_positions)
+
+        # s2 corrections
+        # S2 top and bottom are corrected separately, and cS2 total is the sum of the two
+        # figure out the map name
+        if len(self.s2_map.map_names) > 1:
+            s2_top_map_name = "map_top"
+            s2_bottom_map_name = "map_bottom"
+        else:
+            s2_top_map_name = "map"
+            s2_bottom_map_name = "map"
+
+        for peak_type in ["", "alt_"]:
+            # For electron lifetime corrections to the S2s,
+            # use lifetimes computed using the main S1.
+            el_string = peak_type + "s2_interaction_" if peak_type == "alt_" else peak_type
+            result[f"{peak_type}elife_correction"] = np.exp(events[f'{el_string}drift_time'] / self.elife)
+
+            # S2(x,y) corrections use the observed S2 positions
+            s2_positions = np.vstack([events[f'{peak_type}s2_x'], events[f'{peak_type}s2_y']]).T
+
+            result[f"{peak_type}cs2_top"] = (events[f'{peak_type}s2_area'] * events[f'{peak_type}s2_area_fraction_top']
+                                             * result[f"{peak_type}elife_correction"]
+                                             / self.s2_map(s2_positions, map_name=s2_top_map_name))
+            result[f"{peak_type}cs2_bottom"] = (events[f'{peak_type}s2_area']
+                                                * (1 - events[f'{peak_type}s2_area_fraction_top'])
+                                                * result[f"{peak_type}elife_correction"]
+                                                / self.s2_map(s2_positions, map_name=s2_bottom_map_name))
+            result[f"{peak_type}cs2"] = result[f"{peak_type}cs2_top"] + result[f"{peak_type}cs2_bottom"]
+            result[f"{peak_type}cs2_area_fraction_top"] = result[f"{peak_type}cs2_top"] / result[f"{peak_type}cs2"]
 
-        # For electron lifetime corrections to the S2s,
-        # use lifetimes computed using the main S1.
-        lifetime_corr = np.exp(events['drift_time'] / self.elife)
-        alt_lifetime_corr = (
-            np.exp((events['alt_s2_interaction_drift_time'])
-                   / self.elife))
-
-        # S2(x,y) corrections use the observed S2 positions
-        s2_positions = np.vstack([events['s2_x'], events['s2_y']]).T
-        alt_s2_positions = np.vstack([events['alt_s2_x'], events['alt_s2_y']]).T
-
-        return dict(
-            time=events['time'],
-            endtime=strax.endtime(events),
-
-            cs1=events['s1_area'] / self.s1_map(event_positions),
-            alt_cs1=events['alt_s1_area'] / self.s1_map(event_positions),
-
-            cs2=(events['s2_area'] * lifetime_corr
-                 / self.s2_map(s2_positions)),
-            alt_cs2=(events['alt_s2_area'] * alt_lifetime_corr
-                     / self.s2_map(alt_s2_positions)))
+        return result
 
 
 @export