Added functional elastic test.

	new file:   tests/elastic_test.cxx
modified : tests / elastic_test.sh

tests/elastic_test.cxx

+#include "nlohmann/json.hpp"
+#include <cmath>
+#include <iostream>
+
+#include "ROOT/RDataFrame.hxx"
+#include "ROOT/RVec.hxx"
+
+#include "Math/Vector3D.h"
+#include "Math/Vector4D.h"
+#include "Math/VectorUtil.h"
+#include "TCanvas.h"
+#include "TLatex.h"
+#include "TStyle.h"
+#include "TSystem.h"
+R__LOAD_LIBRARY(libMathMore.so)
+R__LOAD_LIBRARY(libGenVector.so)
+
+R__LOAD_LIBRARY(libfmt.so)
+#include "fmt/core.h"
+
+#include "THStack.h"
+
+#ifdef __cpp_lib_filesystem
+#include <filesystem>
+namespace fs = std::filesystem;
+#else
+#include <experimental/filesystem>
+namespace fs = std::experimental::filesystem;
+#endif
+
+using RDFNode = decltype(ROOT::RDataFrame{0}.Filter(""));
+using Histo1DProxy =
+    decltype(ROOT::RDataFrame{0}.Histo1D(ROOT::RDF::TH1DModel{"", "", 128u, 0., 0.}, ""));
+
+struct RDFInfo {
+  RDFNode&          df;
+  const std::string title;
+  RDFInfo(RDFNode& df, std::string_view title) : df{df}, title{title} {}
+};
+
+constexpr const double M_P     = .938272;
+constexpr const double M_P2    = M_P * M_P;
+constexpr const double M_pion  = 0.139;
+constexpr const double M_pion2 = M_pion * M_pion;
+constexpr const double M_e     = .000511;
+
+// =================================================================================
+// Cuts
+// =================================================================================
+std::string goodTrackSHMS = "P.gtr.dp > -10 && P.gtr.dp < 22";
+std::string goodTrackHMS  = "H.gtr.dp > -8 && H.gtr.dp < 8";
+
+std::string piCutSHMS = "P.cal.etottracknorm<1.0";
+//std::string piCutSHMS = "P.aero.npeSum > 1.0 && P.cal.eprtracknorm < 0.2 && P.cal.etottracknorm<1.0";
+
+std::string eCutHMS = "H.cal.etottracknorm > 0.50 && H.cal.etottracknorm < 2. && "
+                      "H.cer.npeSum > 1.";
+
+std::string epiCut = "P.aero.npeSum > 1.0 && P.cal.eprtracknorm < 0.2 && "
+                     "H.cer.npeSum > 1.0 && H.cal.etottracknorm > 0.6 && "
+                     "H.cal.etottracknorm < 2.0 && P.cal.etottracknorm<1.0";
+
+using Pvec3D = ROOT::Math::XYZVector;
+using Pvec4D = ROOT::Math::PxPyPzMVector;
+
+// =================================================================================
+// reconstruction
+// =================================================================================
+auto p_pion = [](double px, double py, double pz) {
+  return Pvec4D{px * 0.996, py * 0.996, pz * 0.996, M_e};
+};
+auto p_electron = [](double px, double py, double pz) {
+  return Pvec4D{px * 0.994, py * 0.994, pz * 0.994, M_e};
+};
+auto t = [](const double Egamma, Pvec4D& jpsi) {
+  Pvec4D beam{0, 0, Egamma, 0};
+  return (beam - jpsi).M2();
+};
+
+bool root_file_exists(std::string rootfile) {
+  bool found_good_file = false;
+  if (!gSystem->AccessPathName(rootfile.c_str())) {
+    TFile file(rootfile.c_str());
+    if (file.IsZombie()) {
+      std::cout << rootfile << " is a zombie.\n";
+      std::cout
+          << " Did your replay finish?  Check that the it is done before running this script.\n";
+      return false;
+      // return;
+    } else {
+      std::cout << " using : " << rootfile << "\n";
+      return true;
+    }
+  }
+  return false;
+}
+
+void elastic_test(int RunNumber = 6012, int nevents = 50000, int prompt = 0, int update = 0,
+                        int default_count_goal = 10000, int redo_timing = 0) {
+
+  // ===============================================================================================
+  // Initialization
+  // ===============================================================================================
+  using nlohmann::json;
+  json j;
+  {
+    std::ifstream json_input_file("db2/run_list_coin.json");
+    try {
+      json_input_file >> j;
+    } catch (json::parse_error) {
+      std::cerr << "error: json file, db2/run_list.json, is incomplete or has broken syntax.\n";
+      std::quick_exit(-127);
+    }
+  }
+
+  auto runnum_str = std::to_string(RunNumber);
+  if (j.find(runnum_str) == j.end()) {
+    std::cout << "Run " << RunNumber << " not found in db2/run_list_coin.json\n";
+    std::cout << "Check that run number and replay exists. \n";
+    std::cout << "If problem persists please contact Sylvester (217-848-0565)\n";
+  }
+  double P0_shms_setting = j[runnum_str]["spectrometers"]["shms_momentum"].get<double>();
+  double P0_shms         = std::abs(P0_shms_setting);
+
+  bool found_good_file = false;
+
+  std::string rootfile =
+      fmt::format("full_online/coin_replay_production_{}_{}.root", RunNumber, nevents);
+  found_good_file = root_file_exists(rootfile.c_str());
+  if (!found_good_file) {
+    rootfile =
+        fmt::format("ROOTfiles_volatile/coin_replay_production_{}_{}.root", RunNumber, nevents);
+    found_good_file = root_file_exists(rootfile.c_str());
+  }
+  if (!found_good_file) {
+    rootfile = fmt::format("ROOTfiles_csv/coin_replay_production_{}_{}.root", RunNumber, nevents);
+    found_good_file = root_file_exists(rootfile.c_str());
+  }
+  if (!found_good_file) {
+    rootfile = fmt::format("ROOTfiles/coin_replay_production_{}_{}.root", RunNumber, nevents);
+    found_good_file = root_file_exists(rootfile.c_str());
+  }
+  if (!found_good_file) {
+    std::cout << " Error: suitable root file not found\n";
+    return;
+  }
+
+  // ===============================================================================================
+  // Dataframe
+  // ===============================================================================================
+
+  ROOT::EnableImplicitMT(24);
+
+  //---------------------------------------------------------------------------
+  // Detector tree
+  ROOT::RDataFrame d("T", rootfile);
+
+  //// SHMS Scaler tree
+  //ROOT::RDataFrame d_sh("TSP", rootfile);
+  //// int N_scaler_events = *(d_sh.Count());
+
+  auto d_coin = d.Filter("fEvtHdr.fEvtType == 4");
+
+  // Good track cuts
+  auto dHMSGoodTrack  = d_coin.Filter(goodTrackHMS);
+  auto dSHMSGoodTrack = d_coin.Filter(goodTrackSHMS);
+  auto dCOINGoodTrack = dHMSGoodTrack.Filter(goodTrackSHMS)
+                            .Define("p_electron", p_electron, {"H.gtr.px", "H.gtr.py", "H.gtr.pz"})
+                            .Define("p_pion", p_pion, {"P.gtr.px", "P.gtr.py", "P.gtr.pz"});
+  // PID cuts
+  auto dHMSEl  = dHMSGoodTrack.Filter(eCutHMS);
+  auto dSHMSEl = dSHMSGoodTrack.Filter(piCutSHMS);
+  auto dCOINEl = dCOINGoodTrack.Filter(eCutHMS + " && " + piCutSHMS);
+                     //.Filter(
+                     //    [=](double npe, double dp) {
+                     //      double p_track = P0_shms * (100.0 + dp) / 100.0;
+                     //      // no cerenkov cut needed when momentum is below 2.8 GeV/c
+                     //      if (p_track < 2.8) {
+                     //        return true;
+                     //      }
+                     //      return npe > 1.0;
+                     //    },
+                     //    {"P.hgcer.npeSum", "P.gtr.dp"});
+
+  // Timing cuts
+  // Find the timing peak
+  // Find the coin peak
+  double coin_peak_center = 0;
+  if (redo_timing) {
+    auto h_coin_time =
+        dCOINEl.Histo1D({"coin_time", "coin_time", 8000, 0, 1000}, "CTime.ePositronCoinTime_ROC2");
+    h_coin_time->DrawClone();
+    int coin_peak_bin = h_coin_time->GetMaximumBin();
+    coin_peak_center  = h_coin_time->GetBinCenter(coin_peak_bin);
+    std::cout << "COINCIDENCE time peak found at: " << coin_peak_center << std::endl;
+  } else {
+    //coin_peak_center = 43.0; // run 7240-7241
+    coin_peak_center = 23.0; // run 6012
+    std::cout << "COINCIDENCE time peak: using pre-calculated value at: " << coin_peak_center
+              << std::endl;
+    ;
+  }
+  // timing cut lambdas
+  // TODO: evaluate timing cut and offset for random background
+  auto timing_cut = [=](double coin_time) { return std::abs(coin_time - coin_peak_center) < 2.; };
+  // anti-timing set to 5x width of regular
+  auto anti_timing_cut = [=](double coin_time) {
+    return std::abs(coin_time - coin_peak_center - 28.) < 10.;
+  };
+
+  // timing counts
+  auto dHMSElInTime  = dHMSEl.Filter(timing_cut, {"CTime.ePositronCoinTime_ROC2"});
+  auto dHMSElRandom  = dHMSEl.Filter(anti_timing_cut, {"CTime.ePositronCoinTime_ROC2"});
+  auto dSHMSElInTime = dSHMSEl.Filter(timing_cut, {"CTime.ePositronCoinTime_ROC2"});
+  auto dSHMSElRandom = dSHMSEl.Filter(anti_timing_cut, {"CTime.ePositronCoinTime_ROC2"});
+
+  auto dCOINElInTime = dCOINEl.Filter(timing_cut, {"CTime.ePiCoinTime_ROC2"});
+  auto dCOINElRandom = dCOINEl.Filter(anti_timing_cut, {"CTime.ePiCoinTime_ROC2"});
+
+  // Output root file
+  //auto out_file =
+  //    new TFile(fmt::format("monitoring/{}/good_csv_counter.root", RunNumber).c_str(), "UPDATE");
+  //out_file->cd();
+
+  // =========================================================================================
+  // Histograms
+  // =========================================================================================
+  // 2D correlations
+  auto hTheta2DNoCuts = d_coin.Histo2D(
+      {"theta2D", "No cuts;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1, .1, 50, -.1, .1},
+      "P.gtr.th", "H.gtr.th");
+  auto hTheta2DTracking = dCOINGoodTrack.Histo2D(
+      {"theta2D", "Cuts: tracking;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1, .1, 50, -.1, .1},
+      "P.gtr.th", "H.gtr.th");
+  auto hTheta2DPID =
+      dCOINEl.Histo2D({"theta2D", "Cuts: tracking+PID;#theta_{SHMS};#theta_{HMS};#counts", 50, -.1,
+                       .1, 50, -.1, .1},
+                      "P.gtr.th", "H.gtr.th");
+  auto hTheta2DTiming =
+      dCOINElInTime.Histo2D({"theta2D", "Cuts: tracking+PID;#theta_{SHMS};#theta_{HMS};#counts", 50,
+                             -.1, .1, 50, -.1, .1},
+                            "P.gtr.th", "H.gtr.th");
+  // timing
+  auto hCoinTimeNoCuts =
+      d_coin.Histo1D({"coin_time.NoCuts", "No Cuts;coin_time;counts", 8000, 0, 1000},
+                     "CTime.ePositronCoinTime_ROC2");
+  auto hCoinTimeTracking = dCOINGoodTrack.Histo1D(
+      {"coin_time.Tracking", "Cuts: Tracking;coin_time;counts", 8000, 0, 1000},
+      "CTime.ePositronCoinTime_ROC2");
+  auto hCoinTimePID =
+      dCOINEl.Histo1D({"coin_time.PID", "Cuts: Tracking+PID;coin_time;counts", 8000, 0, 1000},
+                      "CTime.ePositronCoinTime_ROC2");
+  auto hCoinTimeTiming = dCOINElInTime.Histo1D(
+      {"coin_time.Timing", "Cuts: Tracking+PID+Timing;coin_time;counts", 8000, 0, 1000},
+      "CTime.ePositronCoinTime_ROC2");
+
+  auto hRandCoinTimePID = dCOINElRandom.Histo1D(
+      {"rand_coin_time.PID", "Cuts: Tracking+PID;coin_time;counts", 8000, 0, 1000},
+      "CTime.ePositronCoinTime_ROC2");
+
+  // P.gtr.dp
+  auto hPdpNoCuts =
+      d_coin.Histo1D({"P.gtr.dp.NoCuts", "No Cuts;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
+  auto hPdpTracking = dSHMSGoodTrack.Histo1D(
+      {"P.gtr.dp.Tracking", "Cuts: Tracking;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
+  auto hPdpPID = dSHMSEl.Histo1D(
+      {"P.gtr.dp.PID", "Cuts: Tracking+PID;#deltap [%];counts", 200, -30, 40}, "P.gtr.dp");
+  auto hPdpTiming = dSHMSElInTime.Histo1D(
+      {"P.gtr.dp.Timing", "Cuts: Tracking+PID+Timing;#deltap [%];counts", 200, -30, 40},
+      "P.gtr.dp");
+  // P.gtr.th
+  auto hPthNoCuts = d_coin.Histo1D(
+      {"P.gtr.th.NoCuts", "No Cuts;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
+  auto hPthTracking = dSHMSGoodTrack.Histo1D(
+      {"P.gtr.th.Tracking", "Cuts: Tracking;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
+  auto hPthPID = dSHMSEl.Histo1D(
+      {"P.gtr.th.PID", "Cuts: Tracking+PID;#theta_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.th");
+  auto hPthTiming = dSHMSElInTime.Histo1D(
+      {"P.gtr.th.Timing", "Cuts: Tracking+PID+Timing;#theta_{SHMS};counts", 200, -0.1, 0.1},
+      "P.gtr.th");
+  // P.gtr.ph
+  auto hPphNoCuts =
+      d_coin.Histo1D({"P.gtr.ph.NoCuts", "No Cuts;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
+  auto hPphTracking = dSHMSGoodTrack.Histo1D(
+      {"P.gtr.ph.Tracking", "Cuts: Tracking;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
+  auto hPphPID = dSHMSEl.Histo1D(
+      {"P.gtr.ph.PID", "Cuts: Tracking+PID;#phi_{SHMS};counts", 200, -0.1, 0.1}, "P.gtr.ph");
+  auto hPphTiming = dSHMSElInTime.Histo1D(
+      {"P.gtr.ph.Timing", "Cuts: Tracking+PID+Timing;#phi_{SHMS};counts", 200, -0.1, 0.1},
+      "P.gtr.ph");
+  // P.gtr.y
+  auto hPyNoCuts =
+      d_coin.Histo1D({"P.gtr.y.NoCuts", "No Cuts;ytar;counts", 200, -10., 10.}, "P.gtr.y");
+  auto hPyTracking = dSHMSGoodTrack.Histo1D(
+      {"P.gtr.y.Tracking", "Cuts: Tracking;ytar;counts", 200, -10., 10.}, "P.gtr.y");
+  auto hPyPID =
+      dSHMSEl.Histo1D({"P.gtr.y.PID", "Cuts: Tracking+PID;ytar;counts", 200, -10., 10.}, "P.gtr.y");
+  auto hPyTiming = dSHMSElInTime.Histo1D(
+      {"P.gtr.y.Timing", "Cuts: Tracking+PID+Timing;ytar;counts", 200, -10., 10.}, "P.gtr.y");
+  // P.cal.etottracknorm
+  auto hPcalEPNoCuts =
+      d_coin.Histo1D({"P.cal.etottracknorm.NoCuts", "No Cuts;SHMS E/P;counts", 200, -.5, 1.5},
+                     "P.cal.etottracknorm");
+  auto hPcalEPTracking = dSHMSGoodTrack.Histo1D(
+      {"P.cal.etottracknorm.Tracking", "Cuts: Tracking;SHMS E/P;counts", 200, -.5, 1.5},
+      "P.cal.etottracknorm");
+  auto hPcalEPPID = dSHMSEl.Histo1D(
+      {"P.cal.etottracknorm.PID", "Cuts: Tracking+PID;SHMS E/P;counts", 200, -.5, 1.5},
+      "P.cal.etottracknorm");
+  auto hPcalEPAll = dCOINElInTime.Histo1D(
+      {"P.cal.etottracknorm.All", "Cuts: Tracking+PID+Coincidence;SHMS E/P;counts", 200, -.5, 1.5},
+      "P.cal.etottracknorm");
+  // P.ngcer.npeSum
+  auto hPcerNpheNoCuts = d_coin.Histo1D(
+      {"P.ngcer.npeSum.NoCuts", "No Cuts;SHMS NGC #phe;counts", 200, -5, 76}, "P.ngcer.npeSum");
+  auto hPcerNpheTracking = dSHMSGoodTrack.Histo1D(
+      {"P.ngcer.npeSum.Tracking", "Cuts: Tracking;SHMS NGC #phe;counts", 200, -5, 76},
+      "P.ngcer.npeSum");
+  auto hPcerNphePID = dSHMSEl.Histo1D(
+      {"P.ngcer.npeSum.PID", "Cuts: Tracking+PID;SHMS NGC #phe;counts", 200, -5, 76},
+      "P.ngcer.npeSum");
+  auto hPcerNpheAll = dCOINElInTime.Histo1D(
+      {"P.ngcer.npeSum.All", "Cuts: Tracking+PID+Coincidence;SHMS NGC #phe;counts", 200, -5, 76},
+      "P.ngcer.npeSum");
+  // P.hgcer.npeSum
+  auto hPhgcerNpheNoCuts = d_coin.Histo1D(
+      {"P.hgcer.npeSum.NoCuts", "No Cuts;SHMS HGC #phe;counts", 200, -5, 76}, "P.hgcer.npeSum");
+  auto hPhgcerNpheTracking = dSHMSGoodTrack.Histo1D(
+      {"P.hgcer.npeSum.Tracking", "Cuts: Tracking;SHMS HGC #phe;counts", 200, -5, 76},
+      "P.hgcer.npeSum");
+  auto hPhgcerNphePID = dSHMSEl.Histo1D(
+      {"P.hgcer.npeSum.PID", "Cuts: Tracking+PID;SHMS HGC #phe;counts", 200, -5, 76},
+      "P.hgcer.npeSum");
+  auto hPhgcerNpheAll = dCOINElInTime.Histo1D(
+      {"P.hgcer.npeSum.All", "Cuts: Tracking+PID+Coincidence;SHMS HGC #phe;counts", 200, -5, 76},
+      "P.hgcer.npeSum");
+  // H.cal.etottracknorm
+  auto hHcalEPNoCuts =
+      d_coin.Histo1D({"H.cal.etottracknorm.NoCuts", "No Cuts;HMS E/P;counts", 200, -.5, 1.5},
+                     "H.cal.etottracknorm");
+  auto hHcalEPTracking = dHMSGoodTrack.Histo1D(
+      {"H.cal.etottracknorm.Tracking", "Cuts: Tracking;HMS E/P;counts", 200, -.5, 1.5},
+      "H.cal.etottracknorm");
+  auto hHcalEPPID = dHMSEl.Histo1D(
+      {"H.cal.etottracknorm.PID", "Cuts: Tracking+PID;HMS E/P;counts", 200, -.5, 1.5},
+      "H.cal.etottracknorm");
+  auto hHcalEPAll = dCOINElInTime.Histo1D(
+      {"H.cal.etottracknorm.All", "Cuts: Tracking+PID+Coincidence;HMS E/P;counts", 200, -.5, 1.5},
+      "H.cal.etottracknorm");
+  // H.cer.npeSum
+  auto hHcerNpheNoCuts = d_coin.Histo1D(
+      {"H.cer.npeSum.NoCuts", "No Cuts;HMS Cer #phe;counts", 200, -1, 15}, "H.cer.npeSum");
+  auto hHcerNpheTracking = dHMSGoodTrack.Histo1D(
+      {"H.cer.npeSum.Tracking", "Cuts: Tracking;HMS Cer #phe;counts", 200, -1, 15}, "H.cer.npeSum");
+  auto hHcerNphePID = dHMSEl.Histo1D(
+      {"H.cer.npeSum.PID", "Cuts: Tracking+PID+Coincidence;HMS Cer #phe;counts", 200, -1, 15},
+      "H.cer.npeSum");
+  auto hHcerNpheAll = dCOINElInTime.Histo1D(
+      {"H.cer.npeSum.PID", "Cuts: Tracking+PID+Coincidence;HMS Cer #phe;counts", 200, -1, 15},
+      "H.cer.npeSum");
+  // H.gtr.dp
+  auto hHdpNoCuts =
+      d_coin.Histo1D({"H.gtr.dp.NoCuts", "No Cuts;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
+  auto hHdpTracking = dHMSGoodTrack.Histo1D(
+      {"H.gtr.dp.Tracking", "Cuts: Tracking;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
+  auto hHdpPID = dHMSEl.Histo1D(
+      {"H.gtr.dp.PID", "Cuts: Tracking+PID;#deltap [%];counts", 200, -30, 40}, "H.gtr.dp");
+  auto hHdpTiming = dHMSElInTime.Histo1D(
+      {"H.gtr.dp.Timing", "Cuts: Tracking+PID+Timing;#deltap [%];counts", 200, -30, 40},
+      "H.gtr.dp");
+  // H.gtr.th
+  auto hHthNoCuts = d_coin.Histo1D(
+      {"H.gtr.th.NoCuts", "No Cuts;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
+  auto hHthTracking = dHMSGoodTrack.Histo1D(
+      {"H.gtr.th.Tracking", "Cuts: Tracking;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
+  auto hHthPID = dHMSEl.Histo1D(
+      {"H.gtr.th.PID", "Cuts: Tracking+PID;#theta_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.th");
+  auto hHthTiming = dHMSElInTime.Histo1D(
+      {"H.gtr.th.Timing", "Cuts: Tracking+PID+Timing;#theta_{HMS};counts", 200, -0.1, 0.1},
+      "H.gtr.th");
+  // H.gtr.ph
+  auto hHphNoCuts =
+      d_coin.Histo1D({"H.gtr.ph.NoCuts", "No Cuts;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
+  auto hHphTracking = dHMSGoodTrack.Histo1D(
+      {"H.gtr.ph.Tracking", "Cuts: Tracking;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
+  auto hHphPID = dHMSEl.Histo1D(
+      {"H.gtr.ph.PID", "Cuts: Tracking+PID;#phi_{HMS};counts", 200, -0.1, 0.1}, "H.gtr.ph");
+  auto hHphTiming = dHMSElInTime.Histo1D(
+      {"H.gtr.ph.Timing", "Cuts: Tracking+PID+Timing;#phi_{HMS};counts", 200, -0.1, 0.1},
+      "H.gtr.ph");
+  // H.gtr.y
+  auto hHyNoCuts =
+      d_coin.Histo1D({"H.gtr.y.NoCuts", "No Cuts;ytar;counts", 200, -10., 10.}, "H.gtr.y");
+  auto hHyTracking = dHMSGoodTrack.Histo1D(
+      {"H.gtr.y.Tracking", "Cuts: Tracking;ytar;counts", 200, -10., 10.}, "H.gtr.y");
+  auto hHyPID =
+      dHMSEl.Histo1D({"H.gtr.y.PID", "Cuts: Tracking+PID;ytar;counts", 200, -10., 10.}, "H.gtr.y");
+  auto hHyTiming = dHMSElInTime.Histo1D(
+      {"H.gtr.y.Timing", "Cuts: Tracking+PID+Timing;ytar;counts", 200, -10., 10.}, "H.gtr.y");
+
+  // scalers
+  //auto total_charge        = d_sh.Max("P.BCM4B.scalerChargeCut");
+  //auto shms_el_real_scaler = d_sh.Max("P.pEL_REAL.scaler");
+  //auto hms_el_real_scaler  = d_sh.Max("P.hEL_REAL.scaler");
+  //auto time_1MHz           = d_sh.Max("P.1MHz.scalerTime");
+  //auto time_1MHz_cut       = d_sh.Max("P.1MHz.scalerTimeCut");
+
+  auto yield_all = d.Count();
+  // 5 timing cut widths worth of random backgrounds
+  auto yield_coin          = d_coin.Count();
+  auto yield_HMSGoodTrack  = dHMSGoodTrack.Count();
+  auto yield_SHMSGoodTrack = dSHMSGoodTrack.Count();
+  auto yield_COINGoodTrack = dCOINGoodTrack.Count();
+  auto yield_HMSEl         = dHMSEl.Count();
+  auto yield_SHMSEl        = dSHMSEl.Count();
+  auto yield_COINEl        = dCOINEl.Count();
+  auto yield_HMSElInTime   = dHMSElInTime.Count();
+  auto yield_HMSElRandom   = dHMSElRandom.Count();
+  auto yield_SHMSElInTime  = dSHMSElInTime.Count();
+  auto yield_SHMSElRandom  = dSHMSElRandom.Count();
+  auto yield_COINElInTime  = dCOINElInTime.Count();
+  auto yield_COINElRandom  = dCOINElRandom.Count();
+  auto yield_coin_raw      = dCOINElInTime.Count();
+  auto yield_coin_random   = dCOINElRandom.Count();
+
+
+  // -------------------------------------
+  // End lazy eval
+  // -------------------------------------
+  auto n_coin_good  = *yield_coin_raw - *yield_coin_random / 5.;
+  auto n_HMSElGood  = *yield_HMSElInTime - *yield_HMSElRandom / 5;
+  auto n_SHMSElGood = *yield_SHMSElInTime - *yield_SHMSElRandom / 5;
+  auto n_COINElGood = *yield_COINElInTime - *yield_COINElRandom / 5;
+
+  hPdpNoCuts->DrawCopy();
+  //std::cout << "  coin COUNTS : " << *(d_coin.Count()) << "\n";
+  //std::cout << " yield_HMSEl : " << *(yield_HMSEl) << "\n";
+  std::cout << " yield_COINEl : " << *(yield_COINEl) << "\n";
+  //std::cout << "  ALL COUNTS : " << *yield_all << "\n";
+  //std::cout << " GOOD COUNTS : " << n_COINElGood << "\n";
+  //
+  if( 4 != (*(yield_COINEl)) ){
+    std::exit(-1);
+  }
+  std::exit(0);
+}

tests/elastic_test.sh

@@ -23,6 +23,4 @@
 singularity help build/Singularity.hcana.simg
 
 
-singularity exec build/Singularity.hcana.simg hcana tests/my_root_script.cxx
-
-echo "woo"
+singularity exec build/Singularity.hcana.simg hcana tests/elastic_test.cxx