diff options
Diffstat (limited to 'src')
-rw-r--r-- | src/Record_Info.h | 1 | ||||
-rw-r--r-- | src/dc.c | 110 | ||||
-rw-r--r-- | src/dc.h | 19 |
3 files changed, 130 insertions, 0 deletions
diff --git a/src/Record_Info.h b/src/Record_Info.h index 91e060f..f3241c6 100644 --- a/src/Record_Info.h +++ b/src/Record_Info.h @@ -214,6 +214,7 @@ typedef struct TriggerInfo { #define TRIG_SPECIAL_RAW 0x00800000 #define TRIG_NCD_MUX 0x01000000 #define TRIG_SOFT_GT 0x02000000 +#define TRIG_MISS 0x04000000 // ------------------------------------------------------------------------------------- // EPED record @@ -25,6 +25,108 @@ #include "zdab_utils.h" #include "zebra.h" #include <gsl/gsl_statistics_double.h> +#include "Record_Info.h" + +#define MAX_PAIRS 10000 + +/* Returns 1 if the event is tagged as failing the FTS cut. The definition for + * the FTS cut comes from the SNOMAN companion: + * + * This tags events where the median time difference of hit PMT pairs + * within 3 m of each other is greater than 6.8 ns. Pairs with large time + * differences (> 25 ns) are not included in the median and at least 15 + * pairs satisfying the distance and time difference criteria are required + * before the event can be tagged. + * + * I used the SNOMAN implementation in flt_fts_cut.for for reference. + * + */ +int is_fts(event *ev) +{ + size_t i, j; + double pmt_dir[3], distance, dt; + double deltat_array[MAX_PAIRS]; + size_t count; + double dtmedian; + + count = 0; + for (i = 0; i < MAX_PMTS; i++) { + /* Only loop over hit PMTs with good calibrations. */ + if (!ev->pmt_hits[i].hit || ev->pmt_hits[i].flags) continue; + + for (j = 0; j < MAX_PMTS; j++) { + /* Only loop over hit PMTs with good calibrations. */ + if (i == j || !ev->pmt_hits[j].hit || ev->pmt_hits[j].flags) continue; + + /* Calculate distance between PMTs. */ + SUB(pmt_dir,pmts[i].pos,pmts[j].pos); + distance = NORM(pmt_dir); + + /* Get the absolute value of the time difference. */ + dt = fabs(ev->pmt_hits[i].t - ev->pmt_hits[j].t); + + if (dt <= 25.0 && distance < 300.0) { + deltat_array[count++] = dt; + if (count >= MAX_PAIRS) goto end; + } + } + } + +end: + + if (count > FTS_COUNT_THRESH) { + gsl_sort(deltat_array,1,count); + dtmedian = gsl_stats_median_from_sorted_data(deltat_array,1,count); + return dtmedian < FTS_MEDIAN_CUT; + } + + return 0; +} + +/* Returns 1 if the event is tagged as an OWL trigger event. The definition for + * the OWL trigger cut comes from the SNOMAN companion: + * + * This tags events where the OWLESUMHI trigger fires. + * + */ +int is_owl_trigger(event *ev) +{ + return (ev->trigger_type & TRIG_OWLE_HI) != 0; +} + +/* Returns 1 if the event is tagged as an OWL event. The definition for the OWL + * cut comes from the SNOMAN companion: + * + * This tags events where the number of OWLs + BUTTs is 3 or more. + * + */ +int is_owl(event *ev) +{ + size_t i; + size_t nhit; + + nhit = 0; + for (i = 0; i < MAX_PMTS; i++) { + if (ev->pmt_hits[i].hit && (pmts[i].pmt_type == PMT_OWL || pmts[i].pmt_type == PMT_BUTT)) + nhit += 1; + } + + return nhit >= 3; +} + +/* Returns 1 if the event fails the ESUM cut. The definition for the ESUM cut + * comes from the SNOMAN companion: + * + * Tags an event if the ESUMLO AND ESUMHI trigger bits are set, AND NONE of + * the following trigger bits are set: (NHIT100LO, NHIT100MED, + * NHIT100HI,NHIT20, NHIT20LB, OWLN, OWLELO, OWLEHI, PULSE_GT, PRESCALE, + * MISS). + * + */ +int is_esum(event *ev) +{ + return (ev->trigger_type & (TRIG_ESUM_LO | TRIG_ESUM_HI | TRIG_NHIT_100_LO | TRIG_NHIT_100_MED | TRIG_NHIT_100_HI | TRIG_NHIT_20 | TRIG_NHIT_20_LB | TRIG_OWLN | TRIG_OWLE_LO | TRIG_OWLE_HI | TRIG_PULSE_GT | TRIG_PRESCALE | TRIG_MISS)) == (TRIG_ESUM_LO | TRIG_ESUM_HI); +} /* Returns the data cleaning bitmask for the event `ev`. */ uint32_t get_dc_word(event *ev, zebraFile *f, zebraBank *bmast, zebraBank *bev) @@ -43,6 +145,14 @@ uint32_t get_dc_word(event *ev, zebraFile *f, zebraBank *bmast, zebraBank *bev) status |= DC_NECK; if (is_flasher(ev)) status |= DC_FLASHER; + if (is_esum(ev)) + status |= DC_ESUM; + if (is_owl(ev)) + status |= DC_OWL; + if (is_owl_trigger(ev)) + status |= DC_OWL_TRIGGER; + if (is_fts(ev)) + status |= DC_FTS; return status; } @@ -28,6 +28,21 @@ #define DC_QVNHIT 0x8 #define DC_NECK 0x10 #define DC_FLASHER 0x20 +#define DC_ESUM 0x40 +#define DC_OWL 0x80 +#define DC_OWL_TRIGGER 0x100 +#define DC_FTS 0x200 + +/* Delta-t threshold for PMT pairs in the FTS cut (ns). */ +#define FTS_DT_THRESH 25.0 +/* Distance threshold for PMT pairs in the FTS cut (cm). */ +#define FTS_DIST_THRESH 300.0 +/* Minimum number of PMT pairs required for the FTS cut. */ +#define FTS_COUNT_THRESH 15 +/* For the FTS cut the median time difference between PMT pairs within 3 meters + * is computed and if this difference is greater than FTS_MEDIAN_CUT the event + * doesn't pass. */ +#define FTS_MEDIAN_CUT 6.8 /* Minimum number of normal PMTs which must be hit to be tagged as an incoming * muon. */ @@ -43,6 +58,10 @@ /* QvNHIT ratio threshold. */ #define QRATIO_THRESHOLD 0.25 +int is_fts(event *ev); +int is_owl_trigger(event *ev); +int is_owl(event *ev); +int is_esum(event *ev); uint32_t get_dc_word(event *ev, zebraFile *f, zebraBank *bmast, zebraBank *bev); int is_muon(event *ev); int junk_cut(zebraFile *f, zebraBank *bmast, zebraBank *bev); |