update is_slot_early() to hopefully catch more flashers

This commit updates is_slot_early() to do the following:

- use median time when checking to see if the time of the PMT hits in the slot
  is early relative to nearby PMTs

- return True if the number of nearby PMTs is less than or equal to the number
  of hit PMTs in the potential flasher slot

- skip potential cross talk hits when finding the nearby hit PMTs by skipping
  hits in the adjacent slots with the same paddle card as the potential flasher

1 files changed, 30 insertions, 15 deletions

diff --git a/src/dc.c b/src/dc.c
index 6195eac..a4dd9f8 100644
--- a/src/dc.c
+++ b/src/dc.c
@@ -649,7 +649,8 @@ int is_breakdown(event *ev)
 
 /* Returns 1 if the average time for the hits in the same card as
  * `flasher_pmt_id` is at least 10 ns earlier than all the hits within 4
- * meters of the potential flasher.
+ * meters of the potential flasher or if there are more hit PMTs in the slot
+ * than within 4 meters.
  *
  * This cut is designed to help flag flashers in the is_flasher() function for
  * events in which the flashing channel may be missing from the event, so there
@@ -657,18 +658,23 @@ int is_breakdown(event *ev)
 int is_slot_early(event *ev, int flasher_pmt_id)
 {
     int i;
-    double t_slot, t_nearby;
+    double t_slot[MAX_PMTS], t_nearby[MAX_PMTS];
+    double t_slot_median, t_nearby_median;
     int n_slot, n_nearby;
     double flasher_pos[3];
     double pmt_dir[3];
     double distance;
+    size_t crate, card, channel;
+    size_t flasher_crate, flasher_card, flasher_channel;
+
+    flasher_crate = flasher_pmt_id/512;
+    flasher_card = (flasher_pmt_id % 512)/32;
+    flasher_channel = flasher_pmt_id % 32;
 
     COPY(flasher_pos,pmts[flasher_pmt_id].pos);
 
     n_slot = 0;
     n_nearby = 0;
-    t_slot = 0.0;
-    t_nearby = 0.0;
     for (i = 0; i < MAX_PMTS; i++) {
         if (!ev->pmt_hits[i].hit || pmts[i].pmt_type != PMT_NORMAL) continue;
 
@@ -721,33 +727,42 @@ int is_slot_early(event *ev, int flasher_pmt_id)
 
         if (i/32 == flasher_pmt_id/32) {
             /* This hit is in the same slot as the potential flasher. */
-            t_slot += ev->pmt_hits[i].ept;
-            n_slot += 1;
+            t_slot[n_slot++] += ev->pmt_hits[i].ept;
             continue;
         }
 
+	crate = i/512;
+	card = (i % 512)/32;
+	channel = i % 32;
+
+        /* Skip PMTs in the same crate that are to the left or right of the
+         * flasher PC since those can often be cross-talk hits. */
+        if (crate == flasher_crate && abs(card - flasher_card) == 1 && (channel/4 == flasher_channel/4)) continue;
+
         /* Calculate the distance from the current channel to the high charge
          * channel. */
         SUB(pmt_dir,pmts[i].pos,flasher_pos);
         distance = NORM(pmt_dir);
 
         if (distance < 400.0) {
-            t_nearby += ev->pmt_hits[i].ept;
-            n_nearby += 1;
+            t_nearby[n_nearby++] += ev->pmt_hits[i].ept;
             continue;
         }
     }
 
-    /* If there are no PMTs within 4 meters of the slot, it's almost certainly
-     * a flasher. */
-    if (n_slot > 4 && n_nearby == 0) return 1;
+    if (n_slot == 0) return 0;
+
+    /* If there are more PMTs in the slot than within 4 meters of the slot,
+     * return it's likely to be a flasher. */
+    if (n_slot >= n_nearby) return 1;
 
-    if (n_slot == 0 || n_nearby == 0) return 0;
+    gsl_sort(t_slot,1,n_slot);
+    t_slot_median = gsl_stats_median_from_sorted_data(t_slot,1,n_slot);
 
-    t_slot /= n_slot;
-    t_nearby /= n_nearby;
+    gsl_sort(t_nearby,1,n_nearby);
+    t_nearby_median = gsl_stats_median_from_sorted_data(t_nearby,1,n_nearby);
 
-    return t_slot < t_nearby - 10.0;
+    return t_slot_median < t_nearby_median - 10.0;
 }
 
 /* Returns 1 if the event is a flasher and 0 if it isn't. The definition of