/* Copyright (c) 2019, Anthony Latorre <tlatorre at uchicago>
 *
 * This program is free software: you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by the Free
 * Software Foundation, either version 3 of the License, or (at your option)
 * any later version.

 * This program is distributed in the hope that it will be useful, but WITHOUT
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
 * more details.

 * You should have received a copy of the GNU General Public License along with
 * this program. If not, see <https://www.gnu.org/licenses/>.
 */

#ifndef DC_H
#define DC_H

#include "event.h"
#include "zebra.h"
#include <stdint.h>

/* Data cleaning bitmasks. */
#define DC_MUON           0x1
#define DC_JUNK           0x2
#define DC_CRATE_ISOTROPY 0x4
#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
#define DC_ITC            0x400
#define DC_BREAKDOWN      0x800

/* Minimum number of hits required for a crate to be considered the source of a
 * breakdown. */
#define MIN_NHIT_BREAKDOWN 256
/* Minimum number of calibrated hits for a crate to calculate the median TAC in
 * the breakdown cut. */
#define MIN_NHIT_CRATE 20
/* Fraction of PMT hits which must have a good TAC value in order to pass the
 * cut. */
#define BREAKDOWN_CAL_FRAC 0.7

/* Length of the sliding window used in the ITC cut (ns). */
#define ITC_TIME_WINDOW 93.0
/* Fraction of hits which must be in a sliding window of 93 ns to pass the ITC
 * cut. */
#define ITC_TIME_FRACTION 0.6

/* Height of PMTs which are considered to be at the "bottom" of the detector
 * for the neck cut. Comes from the SNOOP definition of a neck event in SNO. */
#define NECK_BOTTOM_HEIGHT -425.0
/* Fraction of normal PMT hits which must be below NECK_BOTTOM_HEIGHT. */
#define NECK_BOTTOM_FRACTION 0.5

/* 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. */
#define MUON_MIN_NHIT 150
/* Minimum number of OWL PMTs which must be hit to be tagged as an incoming
 * muon. */
#define MUON_MIN_OWL_NHIT 5
/* Maximum distance between OWL and normal PMTs when computing the median time
 * and charge. */
#define MUON_OWL_NEARBY_DISTANCE 300.0

/* QvNHIT ratio threshold.
 *
 * Note: This is defined as 0.25 in the SNO documentation and in the code, but
 * there is a major bug in the implementation. On line 174 of flt_qnhit_cut.for
 * it does:
 *
 *     * check if pmt terminator is blown
 *     if( .not. btest(ccc_info(KCCC_STATUS,iccc), KCCC_B_75OHM) )then
 *        qhl_work(nhit_cal)=  qhl
 *     else
 *        qhl_work(nhit_cal)=  qhl/2.0
 *     endif
 *
 * I *think* it's always evaluating the qhl/2.0 branch even though the docs and
 * Neil Mccauley's thesis don't say that.
 *
 * I verified this by looking at the QHL over nhit ratio for neutrons in Monte
 * Carlo which peaked around 0.9 whereas the ratio for the N16 events in Neil
 * Mccauley's thesis peaks around 0.5. */
#define QRATIO_THRESHOLD 0.5

int is_itc(event *ev);
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);
int crate_isotropy(event *ev);
int qvnhit(event *ev);
double get_neck_tube_cut_time_diff(event *ev);
int is_neck_event(event *ev);
int is_breakdown(event *ev);
int is_flasher(event *ev);

#endif