Add fiducial & BDT score cut options to EcalPreselectionSkimmer (#1456)

LDMX-Software · Sep 16, 2024 · 5cfe49e · 5cfe49e
1 parent a546306
commit 5cfe49e
Show file tree

Hide file tree

Showing 5 changed files with 58 additions and 27 deletions.
diff --git a/DQM/src/DQM/EcalVetoResults.cxx b/DQM/src/DQM/EcalVetoResults.cxx
@@ -17,7 +17,7 @@ void EcalVetoResults::analyze(const framework::Event &event) {
       event.getObject<ldmx::EcalVetoResult>(ecal_veto_name_, ecal_veto_pass_)};
 
   histograms_.fill("bdt_disc", veto.getDisc());
-  histograms_.fill("bdt_disc_log", -log(1 - veto.getDisc()));
+  histograms_.fill("bdt_disc_log", -std::log10(1 - veto.getDisc()));
   histograms_.fill("fiducial", veto.getFiducial());
 
   return;

diff --git a/Recon/include/Recon/Skims/EcalPreselectionSkimmer.h b/Recon/include/Recon/Skims/EcalPreselectionSkimmer.h
@@ -21,7 +21,7 @@ class EcalPreselectionSkimmer : public framework::Producer {
   EcalPreselectionSkimmer(const std::string &name, framework::Process &process);
 
   /** Destructor */
-  ~EcalPreselectionSkimmer() = default;
+  virtual ~EcalPreselectionSkimmer() = default;
 
   // Configure this processor with a set of parameters passed
   virtual void configure(framework::config::Parameters &) final;
@@ -47,17 +47,26 @@ class EcalPreselectionSkimmer : public framework::Producer {
   /// Max value for num readout hits
   int n_readout_hits_max_;
   /// Max value for shower rms
-  int shower_rms_max_;
+  double shower_rms_max_;
   /// Max value for shower rms in Y
-  int shower_y_std_max_;
+  double shower_y_std_max_;
   /// Max value for shower rms in X
-  int shower_x_std_max_;
+  double shower_x_std_max_;
   /// Max value for maximal cell deposition
   double max_cell_dep_max_;
   /// Max value for std layer hits
   int std_layer_hit_max_;
   /// Max value for num straight tracks
   int n_straight_tracks_max_;
+  /// Min value for the BDT disc variable
+  double bdt_disc_min_;
+  /**
+   * Level of interest in fiducial events
+   *  0: don't care if it's fiducial or not,
+   *  1: keep fiducial events only,
+   *  2: keep non-fid events only
+   */
+  int fiducial_level_;
 
 };  // EcalPreselectionSkimmer
 }  // namespace recon

diff --git a/Recon/python/ecalPreselectionSkimmer.py b/Recon/python/ecalPreselectionSkimmer.py
@@ -17,18 +17,25 @@
     Max value for ecal back energy
  n_readout_hits_max: int
     Max value for num readout hits
- shower_rms_max: int
+ shower_rms_max: double 
     Max value for shower rms
- shower_y_std_max: int
+ shower_y_std_max: double
    Max value for shower rms in Y
-shower_x_std_max: int
+shower_x_std_max: double
     Max value for shower rms in X
 max_cell_dep_max: double
     Max value for maximal cell deposition
 std_layer_hit_max: int
     Max value for std layer hits
 n_straight_tracks_max: int
     Max value for num straight tracks
+bdt_disc_min: double
+    Min value for the BDT disc variable
+fiducial_level: int
+    0: don't care if it's fiducial or not, 
+    1: keep fiducial events only, 
+    2: keep non-fid events only
+
 
 
 Examples
@@ -52,9 +59,12 @@ def __init__(self, name = "ecalPreselectionSkimmer") :
         self.summed_tight_iso_max = 9999.
         self.ecal_back_energy_max = 9999.
         self.n_readout_hits_max = 9999
-        self.shower_rms_max = 9999
-        self.shower_y_std_max = 9999
-        self.shower_x_std_max = 9999
+        self.shower_rms_max = 9999.
+        self.shower_y_std_max = 9999.
+        self.shower_x_std_max = 9999.
         self.max_cell_dep_max = 9999.
         self.std_layer_hit_max = 9999
         self.n_straight_tracks_max = 9999
+        self.bdt_disc_min = 0.
+        self.fiducial_level = 0
+
diff --git a/Recon/src/Recon/Skims/EcalPreselectionSkimmer.cxx b/Recon/src/Recon/Skims/EcalPreselectionSkimmer.cxx
@@ -21,32 +21,42 @@ void EcalPreselectionSkimmer::configure(framework::config::Parameters &ps) {
   ecal_back_energy_max_ =
       ps.getParameter<double>("ecal_back_energy_max");  // MeV
   n_readout_hits_max_ = ps.getParameter<int>("n_readout_hits_max");
-  shower_rms_max_ = ps.getParameter<int>("shower_rms_max");
-  shower_y_std_max_ = ps.getParameter<int>("shower_y_std_max");
-  shower_x_std_max_ = ps.getParameter<int>("shower_x_std_max");
+  shower_rms_max_ = ps.getParameter<double>("shower_rms_max");
+  shower_y_std_max_ = ps.getParameter<double>("shower_y_std_max");
+  shower_x_std_max_ = ps.getParameter<double>("shower_x_std_max");
   max_cell_dep_max_ = ps.getParameter<double>("max_cell_dep_max");  // MeV
   std_layer_hit_max_ = ps.getParameter<int>("std_layer_hit_max");
   n_straight_tracks_max_ = ps.getParameter<int>("n_straight_tracks_max");
+  bdt_disc_min_ = ps.getParameter<double>("bdt_disc_min");
+  fiducial_level_ = ps.getParameter<int>("fiducial_level");
 
   return;
 }
 
 void EcalPreselectionSkimmer::produce(framework::Event &event) {
   bool passedPreselection{false};
+  bool fiducialDecision{true};
   const auto &ecalVeto{
       event.getObject<ldmx::EcalVetoResult>(ecal_veto_name_, ecal_veto_pass_)};
 
+  // Boolean to if we skim for fiducial / nonfiducial
+  fiducialDecision = (fiducial_level_ == 0 ||
+                      (fiducial_level_ == 1 && ecalVeto.getFiducial()) ||
+                      (fiducial_level_ == 2 && !ecalVeto.getFiducial()));
+
   // Boolean to check if we pass preselection
-  passedPreselection = (ecalVeto.getSummedDet() < summed_det_max_) &&
-                       (ecalVeto.getSummedTightIso() < summed_tight_iso_max_) &&
-                       (ecalVeto.getEcalBackEnergy() < ecal_back_energy_max_) &&
-                       (ecalVeto.getNReadoutHits() < n_readout_hits_max_) &&
-                       (ecalVeto.getShowerRMS() < shower_rms_max_) &&
-                       (ecalVeto.getYStd() < shower_y_std_max_) &&
-                       (ecalVeto.getXStd() < shower_x_std_max_) &&
-                       (ecalVeto.getMaxCellDep() < max_cell_dep_max_) &&
-                       (ecalVeto.getStdLayerHit() < std_layer_hit_max_) &&
-                       (ecalVeto.getNStraightTracks() < n_straight_tracks_max_);
+  passedPreselection =
+      (ecalVeto.getSummedDet() < summed_det_max_) &&
+      (ecalVeto.getSummedTightIso() < summed_tight_iso_max_) &&
+      (ecalVeto.getEcalBackEnergy() < ecal_back_energy_max_) &&
+      (ecalVeto.getNReadoutHits() < n_readout_hits_max_) &&
+      (ecalVeto.getShowerRMS() < shower_rms_max_) &&
+      (ecalVeto.getYStd() < shower_y_std_max_) &&
+      (ecalVeto.getXStd() < shower_x_std_max_) &&
+      (ecalVeto.getMaxCellDep() < max_cell_dep_max_) &&
+      (ecalVeto.getStdLayerHit() < std_layer_hit_max_) &&
+      (ecalVeto.getNStraightTracks() < n_straight_tracks_max_) &&
+      (ecalVeto.getDisc() > bdt_disc_min_) && fiducialDecision;
 
   // Tell the skimmer to keep or drop the event based on whether preselection
   // passed

diff --git a/Recon/test/ecal_preselection_skim.py b/Recon/test/ecal_preselection_skim.py
@@ -21,12 +21,14 @@
 ecal_pres_skimmer.summed_tight_iso_max = 9999. 
 ecal_pres_skimmer.ecal_back_energy_max = 9999. 
 ecal_pres_skimmer.n_readout_hits_max = 9999
-ecal_pres_skimmer.shower_rms_max = 9999 
-ecal_pres_skimmer.shower_y_std_max = 9999 
-ecal_pres_skimmer.shower_x_std_max = 9999 
+ecal_pres_skimmer.shower_rms_max = 9999. 
+ecal_pres_skimmer.shower_y_std_max = 9999. 
+ecal_pres_skimmer.shower_x_std_max = 9999.
 ecal_pres_skimmer.max_cell_dep_max = 9999. 
 ecal_pres_skimmer.std_layer_hit_max = 9999 
 ecal_pres_skimmer.n_straight_tracks_max = 9999
+ecal_pres_skimmer.bdt_disc_min = 0.
+ecal_pres_skimmer.fiducial_level = 0
 '''
 
 p.sequence =[ecal_pres_skimmer]