diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_common_functions.h b/benchmarks/barrel_ecal/scripts/emcal_barrel_common_functions.h index b1359567f8e2434cb2880f675c5ee1086d98858a..34e397abfbbdcb570ae6d783099bb7f00a6e7a47 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_common_functions.h +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_common_functions.h @@ -3,6 +3,10 @@ // M. Scott 05/2021 ////////////////////////// +#include <cctype> +#include "TH1.h" + +// Returns particle pdgID and mass in [GeV] std::tuple <int, double> extract_particle_parameters(std::string particle_name) { if (particle_name == "electron") return std::make_tuple(11, 0.51099895e-3); if (particle_name == "photon") return std::make_tuple(22, 0.0); @@ -15,4 +19,18 @@ std::tuple <int, double> extract_particle_parameters(std::string particle_name) std::cout << "wrong particle name" << std::endl; abort(); +} + +// Returns Environment Variables +std::string getEnvVar( std::string const & key ){ + char * val = getenv( key.c_str() ); + return val == NULL ? std::string("") : std::string(val); +} +// Added detetcor name to title +void addDetectorName(std::string name, TH1 * inhist){ + std::string newName = inhist -> GetTitle(); + for (auto& x : name){ + x = toupper(x); + } + inhist -> SetTitle((name + " " + newName).c_str()); } \ No newline at end of file diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_electrons_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_electrons_analysis.cxx index 6b52660ae1d9ac0d952fe1eccc2cf07812098fbf..22cb118ce62847bd6ebef8e313537edbd04e3f80 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_electrons_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_electrons_analysis.cxx @@ -16,11 +16,15 @@ #include "TF1.h" #include "TH1D.h" +#include "emcal_barrel_common_functions.h" +#include <cctype> + using ROOT::RDataFrame; using namespace ROOT::VecOps; void emcal_barrel_electrons_analysis(const char* input_fname = "sim_output/sim_emcal_barrel_uniform_electrons.root") { + input_fname = "temp_pions_electrons.root"; // Setting for graphs gROOT->SetStyle("Plain"); gStyle->SetOptFit(1); @@ -35,6 +39,16 @@ void emcal_barrel_electrons_analysis(const char* input_fname = "sim_output/sim_e ROOT::EnableImplicitMT(); ROOT::RDataFrame d0("events", input_fname); + // Environment Variables + std::string detector_path = ""; + std::string detector_name = "topside"; + if(std::getenv("DETECTOR_PATH")) { + detector_path = std::getenv("DETECTOR_PATH"); + } + if(std::getenv("JUGGLER_DETECTOR")) { + detector_name = std::getenv("JUGGLER_DETECTOR"); + } + // Thrown Energy [GeV] auto Ethr = [](std::vector<dd4pod::Geant4ParticleData> const& input) { auto p = input[2]; @@ -79,6 +93,10 @@ void emcal_barrel_electrons_analysis(const char* input_fname = "sim_output/sim_e {"hfsam", "Sampling Fraction; Sampling Fraction; Events", 100, 0.0, 0.1}, "fsam"); + addDetectorName(detector_name, hEthr.GetPtr()); + addDetectorName(detector_name, hEsim.GetPtr()); + addDetectorName(detector_name, hfsam.GetPtr()); + // Event Counts auto nevents_thrown = d1.Count(); std::cout << "Number of Thrown Events: " << (*nevents_thrown) << "\n"; diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx index 33a31897fe97828d1acc0b65a2cf5f177ced9f7f..b0184884c30349fdd1fcb196e927079e84cced2b 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_particles_analysis.cxx @@ -51,6 +51,16 @@ void emcal_barrel_particles_analysis(std::string particle_name = "electron") std::string input_fname = fmt::format("sim_output/sim_emcal_barrel_{}.root", particle_name); ROOT::RDataFrame d0("events", input_fname); + // Environment Variables + std::string detector_path = ""; + std::string detector_name = "topside"; + if(std::getenv("DETECTOR_PATH")) { + detector_path = std::getenv("DETECTOR_PATH"); + } + if(std::getenv("JUGGLER_DETECTOR")) { + detector_name = std::getenv("JUGGLER_DETECTOR"); + } + // Thrown Energy [GeV] auto Ethr = [](std::vector<dd4pod::Geant4ParticleData> const& input) { auto p = input[2]; @@ -94,6 +104,10 @@ void emcal_barrel_particles_analysis(std::string particle_name = "electron") {"hfsam", "Sampling Fraction; Sampling Fraction; Events", 200, 0.0, 0.05}, "fsam"); + addDetectorName(detector_name, hEthr.GetPtr()); + addDetectorName(detector_name, hEsim.GetPtr()); + addDetectorName(detector_name, hfsam.GetPtr()); + // Event Counts auto nevents_thrown = d1.Count(); std::cout << "Number of Thrown Events: " << (*nevents_thrown) << "\n"; diff --git a/benchmarks/barrel_ecal/scripts/emcal_barrel_pions_electrons_analysis.cxx b/benchmarks/barrel_ecal/scripts/emcal_barrel_pions_electrons_analysis.cxx index 70b1043092cc1ed7608b7f8b35aae8e5662c9440..3bf036cde29b3900cee9465c9e51b2cd178046ae 100644 --- a/benchmarks/barrel_ecal/scripts/emcal_barrel_pions_electrons_analysis.cxx +++ b/benchmarks/barrel_ecal/scripts/emcal_barrel_pions_electrons_analysis.cxx @@ -24,7 +24,7 @@ R__LOAD_LIBRARY(libfmt.so) #include "DD4hep/Detector.h" #include "DDG4/Geant4Data.h" #include "DDRec/CellIDPositionConverter.h" - +#include "emcal_barrel_common_functions.h" using ROOT::RDataFrame; using namespace ROOT::VecOps; @@ -45,7 +45,7 @@ void emcal_barrel_pions_electrons_analysis(const char* input_fname1 = "sim_outpu ROOT::EnableImplicitMT(); ROOT::RDataFrame d0("events", {input_fname1, input_fname2}); - //Using the detector layers + // Environment Variables std::string detector_path = ""; std::string detector_name = "topside"; if(std::getenv("DETECTOR_PATH")) { @@ -55,6 +55,7 @@ void emcal_barrel_pions_electrons_analysis(const char* input_fname1 = "sim_outpu detector_name = std::getenv("JUGGLER_DETECTOR"); } + // Using the detector layers dd4hep::Detector& detector = dd4hep::Detector::getInstance(); detector.fromCompact(fmt::format("{}/{}.xml", detector_path,detector_name)); //dd4hep::rec::CellIDPositionConverter cellid_converter(detector); @@ -156,48 +157,66 @@ void emcal_barrel_pions_electrons_analysis(const char* input_fname1 = "sim_outpu auto d_pim_rej = d_pim.Filter("Esim_front > " + to_string(rejectCut)); // Define Histograms - auto hEthr = d1.Histo1D({"hEthr", "Thrown Energy; Thrown Energy [GeV]; Events", 100, 0.0, 7.5}, "Ethr"); - auto hNhits = d1.Histo1D({"hNhits", "Number of hits per events; Number of hits; Events", 100, 0.0, 2000.0}, "nhits"); - auto hEsim = d1.Histo1D({"hEsim", "Energy Deposit; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim"); - auto hfsam = d1.Histo1D({"hfsam", "Sampling Fraction; Sampling Fraction; Events", 100, 0.0, 0.1}, "fsam"); + auto hEthr = d1.Histo1D({"hEthr", "Thrown Energy; Thrown Energy [GeV]; Events", 100, 0.0, 7.5}, "Ethr"); + auto hNhits = d1.Histo1D({"hNhits", "Number of hits per events; Number of hits; Events", 100, 0.0, 2000.0}, "nhits"); + auto hEsim = d1.Histo1D({"hEsim", "Energy Deposit; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim"); + auto hfsam = d1.Histo1D({"hfsam", "Sampling Fraction; Sampling Fraction; Events", 100, 0.0, 0.1}, "fsam"); auto hEsim_front = d1.Histo1D({"hEsim_front", "Energy Deposit Front; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim_front"); + addDetectorName(detector_name, hEthr.GetPtr()); + addDetectorName(detector_name, hNhits.GetPtr()); + addDetectorName(detector_name, hfsam.GetPtr()); + addDetectorName(detector_name, hEsim_front.GetPtr()); auto hEsim_ele = d_ele.Histo1D({"hEsim_ele", "Energy Deposit Electron; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim"); auto hEsim_ele_front = d_ele.Histo1D({"hEsim_ele_front", "Energy Deposit Front Electron; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim_front"); auto hEsim_pim = d_pim.Histo1D({"hEsim_pim", "Energy Deposit Electron; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim"); auto hEsim_pim_front = d_pim.Histo1D({"hEsim_pim_front", "Energy Deposit Front Pion-; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim_front"); + addDetectorName(detector_name, hEsim_ele.GetPtr()); + addDetectorName(detector_name, hEsim_pim_front.GetPtr()); + addDetectorName(detector_name, hEsim_pim.GetPtr()); + addDetectorName(detector_name, hEsim_pim_front.GetPtr()); auto hEsim_ele_front_rej = d_ele_rej.Histo1D({"hEsim_ele_front_rej", "Energy Deposit Front Electron; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim_front"); auto hEsim_pim_front_rej = d_pim_rej.Histo1D({"hEsim_pim_front_rej", "Energy Deposit Front Pion-; Energy Deposit [GeV]; Events", 10, 0.0, 0.05}, "Esim_front"); + addDetectorName(detector_name, hEsim_ele_front_rej.GetPtr()); + addDetectorName(detector_name, hEsim_pim_front_rej.GetPtr()); auto hEpvp_ele = d_ele.Histo2D({"hEpvp_ele", "Energy Deposit 1st 4 Layers/P vs P; E/P; P [GeV]", 20, 0.0, 0.01, 20, 1.0, 7.0}, "EOverP", "Pthr"); auto hEpvp_pim = d_pim.Histo2D({"hEpvp_pim", "Energy Deposit 1st 4 Layers/P vs P; E/P; P [GeV]", 20, 0.0, 0.01, 20, 1.0, 7.0}, "EOverP", "Pthr"); + addDetectorName(detector_name, hEpvp_ele.GetPtr()); + addDetectorName(detector_name, hEpvp_pim.GetPtr()); TH1D* hElePurity_initial = (TH1D *)hEsim_ele -> Clone(); hElePurity_initial -> Divide(hEsim.GetPtr()); hElePurity_initial -> SetTitle("Electron/Pion Rejection"); + addDetectorName(detector_name, hElePurity_initial); TH1D* hElePurity_ele = (TH1D *)hEsim_ele_front -> Clone(); hElePurity_ele -> Divide(hEsim_front.GetPtr()); hElePurity_ele -> SetTitle("Electron/Pion Rejection : Electron"); + addDetectorName(detector_name, hElePurity_ele); TH1D* hElePurity_pim = (TH1D *)hEsim_pim_front -> Clone(); hElePurity_pim -> Divide(hEsim_front.GetPtr()); hElePurity_pim -> SetTitle("Electron/Pion Rejection : Pion Minus"); + addDetectorName(detector_name, hElePurity_pim); // Rejection TH1D* hElePurity_rej = (TH1D *)hEsim_ele_front_rej -> Clone(); hElePurity_rej -> Add(hEsim_pim_front_rej.GetPtr()); hElePurity_rej -> SetTitle("Electron/Pion Rejection"); + addDetectorName(detector_name, hElePurity_rej); TH1D* hElePurity_ele_rej = (TH1D *)hEsim_ele_front_rej -> Clone(); hElePurity_ele_rej -> Divide(hElePurity_rej); hElePurity_ele_rej -> SetTitle("Electron/Pion Rejection : Electron"); + addDetectorName(detector_name, hElePurity_ele_rej); TH1D* hElePurity_pim_rej = (TH1D *)hEsim_pim_front_rej -> Clone(); hElePurity_pim_rej -> Divide(hElePurity_rej); hElePurity_pim_rej -> SetTitle("Electron/Pion Rejection : Pi Minus"); + addDetectorName(detector_name, hElePurity_pim_rej); // Event Counts auto nevents_thrown = d1.Count();