From 38c3b530ef2e356d9bab8a57abaaeae135e01e92 Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Fri, 6 Oct 2023 10:15:31 +0200 Subject: [PATCH 1/7] Added methods for SimCaloHits --- .../dataframe/FCCAnalyses/CaloNtupleizer.h | 16 +++ analyzers/dataframe/src/CaloNtupleizer.cc | 98 +++++++++++++++++++ 2 files changed, 114 insertions(+) diff --git a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h index 1b9b9a8291..a0b310f0f7 100644 --- a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h +++ b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h @@ -7,6 +7,7 @@ #include "ROOT/RVec.hxx" #include "edm4hep/CalorimeterHitData.h" +#include "edm4hep/SimCalorimeterHitData.h" #include "edm4hep/ClusterData.h" #include "edm4hep/MCParticleData.h" @@ -31,6 +32,21 @@ struct sel_layers { }; +// SIM calo hits (single cells) +ROOT::VecOps::RVec getSimCaloHit_x (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_y (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_z (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_phi (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_phiBin (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_theta (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_etaBin (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps::RVec& in); + + + // calo hits (single cells) ROOT::VecOps::RVec getCaloHit_x (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getCaloHit_y (const ROOT::VecOps::RVec& in); diff --git a/analyzers/dataframe/src/CaloNtupleizer.cc b/analyzers/dataframe/src/CaloNtupleizer.cc index c162cf112b..3589a900dd 100644 --- a/analyzers/dataframe/src/CaloNtupleizer.cc +++ b/analyzers/dataframe/src/CaloNtupleizer.cc @@ -33,7 +33,105 @@ ROOT::VecOps::RVec sel_layers::operator() (const R } return res; } +// SIM calo hit +ROOT::VecOps::RVec getSimCaloHit_x (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.position.x); + } + return result; +} + +ROOT::VecOps::RVec getSimCaloHit_y (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.position.y); + } + return result; +} + +ROOT::VecOps::RVec getSimCaloHit_z (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.position.z); + } + return result; +} + +ROOT::VecOps::RVec getSimCaloHit_phi (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + TVector3 t3; + t3.SetXYZ(p.position.x, p.position.y, p.position.z); + result.push_back(t3.Phi()); + } + return result; +} + +//ROOT::VecOps::RVec getSimCaloHit_phiBin (const ROOT::VecOps::RVec& in){ +// ROOT::VecOps::RVec result; +// for (auto & p: in){ +// dd4hep::DDSegmentation::CellID cellId = p.cellID; +// result.push_back(m_decoder->getSim(cellId, "phi")); +// } +// return result; +//} + +ROOT::VecOps::RVec getSimCaloHit_theta (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + TVector3 t3; + t3.SetXYZ(p.position.x, p.position.y, p.position.z); + result.push_back(t3.Theta()); + } + return result; +} + +ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + TVector3 t3; + t3.SetXYZ(p.position.x, p.position.y, p.position.z); + result.push_back(t3.Eta()); + } + return result; +} + +//ROOT::VecOps::RVec getSimCaloHit_etaBin (const ROOT::VecOps::RVec& in){ +// ROOT::VecOps::RVec result; +// for (auto & p: in){ +// dd4hep::DDSegmentation::CellID cellId = p.cellID; +// result.push_back(m_decoder->get(cellId, "eta")); +// } +// return result; +//} + +ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.energy); + } + return result; +} + +ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + dd4hep::DDSegmentation::CellID cellId = p.cellID; + result.push_back(m_decoder->get(cellId, "layer")); + } + return result; +} +ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + TVector3 t3; + t3.SetXYZ(p.position.x, p.position.y, p.position.z); + result.push_back(t3); + } + return result; +} // calo hit ROOT::VecOps::RVec getCaloHit_x (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; From 726f776419bf257f7054b46bf0f5d28449f9edd5 Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Thu, 12 Oct 2023 10:48:38 +0200 Subject: [PATCH 2/7] added scepcal folder --- SCEPCal_plots/ntuplizer.py | 81 +++++++++++++++++++ .../dataframe/FCCAnalyses/CaloNtupleizer.h | 2 + analyzers/dataframe/src/CaloNtupleizer.cc | 19 +++++ 3 files changed, 102 insertions(+) create mode 100644 SCEPCal_plots/ntuplizer.py diff --git a/SCEPCal_plots/ntuplizer.py b/SCEPCal_plots/ntuplizer.py new file mode 100644 index 0000000000..77e5a3407a --- /dev/null +++ b/SCEPCal_plots/ntuplizer.py @@ -0,0 +1,81 @@ +#!/usr/bin/env python3 +import ROOT +import argparse + +ROOT.EnableImplicitMT() +print ("Load cxx analyzers ... ") +ROOT.gSystem.Load("libFCCAnalyses") + +ROOT.gErrorIgnoreLevel = ROOT.kFatal +_fcc = ROOT.dummyLoader +ROOT.gInterpreter.Declare("using namespace FCCAnalyses;") + +#to add layer info - not working +#ROOT.CaloNtupleizer.loadGeometry('/afs/cern.ch/work/f/fcetorel/private/work2/SCEPCAL/Detector/SCEPCAL/compact/SCEPCAL.xml', 'SCEPCALreadout') + +#/afs/cern.ch/work/f/fcetorel/private/work2/SCEPCAL/SCEPCALsim/SCEPCALsimG4Components/test/data/ + +parser = argparse.ArgumentParser() +parser.add_argument("-inputFilesName", help = "name of the input rootfiles", type = str) +parser.add_argument("-baseFolder",help = "folder for the rootfiles", type = str) +args = parser.parse_args() + + + +print('Create RDataFrame ...') +df = ROOT.RDataFrame('events', args.baseFolder + "/" + args.inputFilesName) +print('Apply selectors and define new branches ...') +df = (df + .Define('SimCaloHit_cellID', 'CaloNtupleizer::getSimCellID (SimCalorimeterHits)') + .Define('SimCaloHit_energy', 'CaloNtupleizer::getSimCaloHit_energy (SimCalorimeterHits)') + .Define('SimCaloHit_r', 'CaloNtupleizer::getSimCaloHit_r(SimCalorimeterHits)') + .Define('SimCaloHit_x', 'CaloNtupleizer::getSimCaloHit_x (SimCalorimeterHits)') + .Define('SimCaloHit_y', 'CaloNtupleizer::getSimCaloHit_y (SimCalorimeterHits)') + .Define('SimCaloHit_z', 'CaloNtupleizer::getSimCaloHit_z (SimCalorimeterHits)') + .Define('SimCaloHit_eta', 'CaloNtupleizer::getSimCaloHit_eta (SimCalorimeterHits)') + .Define('SimCaloHit_phi', 'CaloNtupleizer::getSimCaloHit_phi (SimCalorimeterHits)') + + + .Define('SimCaloHitC_cellID', 'CaloNtupleizer::getSimCellID (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_energy', 'CaloNtupleizer::getSimCaloHit_energy (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_r', 'CaloNtupleizer::getSimCaloHit_r (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_x', 'CaloNtupleizer::getSimCaloHit_x (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_y', 'CaloNtupleizer::getSimCaloHit_y (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_z', 'CaloNtupleizer::getSimCaloHit_z (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_eta', 'CaloNtupleizer::getSimCaloHit_eta (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_phi', 'CaloNtupleizer::getSimCaloHit_phi (SimCalorimeterHitsCherenkov)') + + .Define('MCparticle_px', 'MCParticle::get_px (GenParticles)') + .Define('MCparticle_py', 'MCParticle::get_py (GenParticles)') + .Define('MCparticle_pz', 'MCParticle::get_pz (GenParticles)') + + ) + +outfilename = 'flatNtupla_'+ args.inputFilesName +print(f'Writing snapshot to disk ... \t{outfilename}') + +df.Snapshot('events', outfilename, + [ + 'SimCaloHit_cellID', + 'SimCaloHit_energy', + 'SimCaloHit_r', + 'SimCaloHit_x', + 'SimCaloHit_y', + 'SimCaloHit_z', + 'SimCaloHit_eta', + 'SimCaloHit_phi', + + 'SimCaloHitC_cellID', + 'SimCaloHitC_energy', + 'SimCaloHitC_r', + 'SimCaloHitC_x', + 'SimCaloHitC_y', + 'SimCaloHitC_z', + 'SimCaloHitC_eta', + 'SimCaloHitC_phi', + + 'MCparticle_px', + 'MCparticle_py', + 'MCparticle_pz', + ] + ) diff --git a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h index a0b310f0f7..8acdc540d2 100644 --- a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h +++ b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h @@ -33,6 +33,8 @@ struct sel_layers { // SIM calo hits (single cells) +ROOT::VecOps::RVec getSimCellID (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_r (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_x (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_y (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_z (const ROOT::VecOps::RVec& in); diff --git a/analyzers/dataframe/src/CaloNtupleizer.cc b/analyzers/dataframe/src/CaloNtupleizer.cc index 3589a900dd..b20b069d2b 100644 --- a/analyzers/dataframe/src/CaloNtupleizer.cc +++ b/analyzers/dataframe/src/CaloNtupleizer.cc @@ -34,6 +34,14 @@ ROOT::VecOps::RVec sel_layers::operator() (const R return res; } // SIM calo hit +ROOT::VecOps::RVec getSimCaloHit_r (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(std::sqrt(p.position.x * p.position.x + p.position.y * p.position.y)); + } + return result; +} + ROOT::VecOps::RVec getSimCaloHit_x (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -106,6 +114,16 @@ ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec getSimCellID (const ROOT::VecOps::RVec& in){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.cellID); + } + return result; +} + + + ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -114,6 +132,7 @@ ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ From f1f00053c58b1d990065b8001a628b2251c6dce4 Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Tue, 7 Nov 2023 11:24:08 +0100 Subject: [PATCH 3/7] add layer info for crystals (very harcoded) --- SCEPCal_plots/ntuplizer.py | 9 +++++++-- analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h | 1 + analyzers/dataframe/src/CaloNtupleizer.cc | 7 +++++++ 3 files changed, 15 insertions(+), 2 deletions(-) diff --git a/SCEPCal_plots/ntuplizer.py b/SCEPCal_plots/ntuplizer.py index 77e5a3407a..f280f9a0e8 100644 --- a/SCEPCal_plots/ntuplizer.py +++ b/SCEPCal_plots/ntuplizer.py @@ -19,14 +19,14 @@ parser.add_argument("-inputFilesName", help = "name of the input rootfiles", type = str) parser.add_argument("-baseFolder",help = "folder for the rootfiles", type = str) args = parser.parse_args() - - +NBITS = 27 print('Create RDataFrame ...') df = ROOT.RDataFrame('events', args.baseFolder + "/" + args.inputFilesName) print('Apply selectors and define new branches ...') df = (df .Define('SimCaloHit_cellID', 'CaloNtupleizer::getSimCellID (SimCalorimeterHits)') + .Define('SimCaloHit_depth', 'CaloNtupleizer::getSimCaloHit_depth(SimCalorimeterHits, %i )'%NBITS) .Define('SimCaloHit_energy', 'CaloNtupleizer::getSimCaloHit_energy (SimCalorimeterHits)') .Define('SimCaloHit_r', 'CaloNtupleizer::getSimCaloHit_r(SimCalorimeterHits)') .Define('SimCaloHit_x', 'CaloNtupleizer::getSimCaloHit_x (SimCalorimeterHits)') @@ -37,6 +37,7 @@ .Define('SimCaloHitC_cellID', 'CaloNtupleizer::getSimCellID (SimCalorimeterHitsCherenkov)') + .Define('SimCaloHitC_depth', 'CaloNtupleizer::getSimCaloHit_depth(SimCalorimeterHits, %i)'%NBITS) .Define('SimCaloHitC_energy', 'CaloNtupleizer::getSimCaloHit_energy (SimCalorimeterHitsCherenkov)') .Define('SimCaloHitC_r', 'CaloNtupleizer::getSimCaloHit_r (SimCalorimeterHitsCherenkov)') .Define('SimCaloHitC_x', 'CaloNtupleizer::getSimCaloHit_x (SimCalorimeterHitsCherenkov)') @@ -45,6 +46,7 @@ .Define('SimCaloHitC_eta', 'CaloNtupleizer::getSimCaloHit_eta (SimCalorimeterHitsCherenkov)') .Define('SimCaloHitC_phi', 'CaloNtupleizer::getSimCaloHit_phi (SimCalorimeterHitsCherenkov)') + .Define('MCparticle_energy', 'MCParticle::get_e (GenParticles)') .Define('MCparticle_px', 'MCParticle::get_px (GenParticles)') .Define('MCparticle_py', 'MCParticle::get_py (GenParticles)') .Define('MCparticle_pz', 'MCParticle::get_pz (GenParticles)') @@ -57,6 +59,7 @@ df.Snapshot('events', outfilename, [ 'SimCaloHit_cellID', + 'SimCaloHit_depth', 'SimCaloHit_energy', 'SimCaloHit_r', 'SimCaloHit_x', @@ -66,6 +69,7 @@ 'SimCaloHit_phi', 'SimCaloHitC_cellID', + 'SimCaloHitC_depth', 'SimCaloHitC_energy', 'SimCaloHitC_r', 'SimCaloHitC_x', @@ -74,6 +78,7 @@ 'SimCaloHitC_eta', 'SimCaloHitC_phi', + 'MCparticle_energy', 'MCparticle_px', 'MCparticle_py', 'MCparticle_pz', diff --git a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h index 8acdc540d2..909b7ac70e 100644 --- a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h +++ b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h @@ -44,6 +44,7 @@ ROOT::VecOps::RVec getSimCaloHit_theta (const ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_etaBin (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in); +ROOT::VecOps::RVec getSimCaloHit_depth (const ROOT::VecOps::RVec& in,const int decodingVal); ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps::RVec& in); diff --git a/analyzers/dataframe/src/CaloNtupleizer.cc b/analyzers/dataframe/src/CaloNtupleizer.cc index b20b069d2b..14d5ad0769 100644 --- a/analyzers/dataframe/src/CaloNtupleizer.cc +++ b/analyzers/dataframe/src/CaloNtupleizer.cc @@ -132,6 +132,13 @@ ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec getSimCaloHit_depth (const ROOT::VecOps::RVec& in,const int decodingVal){ + ROOT::VecOps::RVec result; + for (auto & p: in){ + result.push_back(p.cellID >> decodingVal & (8-1) ); + } + return result; +} ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; From 31894b07059df38f2844a7a9e5fa735aaf3c21e0 Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Tue, 7 Nov 2023 16:40:11 +0100 Subject: [PATCH 4/7] added plotter --- SCEPCal_plots/simplePlotter.py | 252 +++++++++++++++++++++++++++++++++ 1 file changed, 252 insertions(+) create mode 100755 SCEPCal_plots/simplePlotter.py diff --git a/SCEPCal_plots/simplePlotter.py b/SCEPCal_plots/simplePlotter.py new file mode 100755 index 0000000000..02168a2723 --- /dev/null +++ b/SCEPCal_plots/simplePlotter.py @@ -0,0 +1,252 @@ +#!/usr/bin/env python +### Produce very simple preliminary plots, no cuts applied, no matching of Cherenkov - Scintillation +import sys, os +import ROOT +from array import array +import argparse + +ROOT.gStyle.SetOptStat(0) +ROOT.gStyle.SetOptFit(1) +ROOT.gStyle.SetOptTitle(0) +ROOT.gStyle.SetLabelSize(0.05,'X') +ROOT.gStyle.SetLabelSize(0.05,'Y') +ROOT.gStyle.SetTitleSize(0.06,'X') +ROOT.gStyle.SetTitleSize(0.06,'Y') +ROOT.gStyle.SetTitleOffset(0.8,'X') +ROOT.gStyle.SetTitleOffset(0.8,'Y') +ROOT.gStyle.SetLegendFont(42) +ROOT.gStyle.SetLegendTextSize(0.038) +ROOT.gStyle.SetPadTopMargin(0.07) +ROOT.gROOT.SetBatch(True) +ROOT.gErrorIgnoreLevel = ROOT.kWarning + +parser = argparse.ArgumentParser(description='Module characterization summary plots') +parser.add_argument("--inputFile", required=True, type=str, help="path to inputfile.root [output of ntuplizer]") +parser.add_argument("--outFolder", required=True, type=str, help="out folder for plots") +args = parser.parse_args() + +# input file should be flat ntupla from ntuplizer +#outdir = "/eos/user/f/fcetorel/www/SCEPCal_Sim/test_plots/newCommitsChecks/" +#inputfile = "flatNtupla_scepcal_100evs_newPull.root" +outdir = args.outFolder +inputfile = args.inputFile +os.makedirs(outdir, exist_ok=True) +os.system("cp /eos/user/f/fcetorel/www/index.php %s"%outdir) + + + +#Open the file - add some branches - create the canvas +d = ROOT.RDataFrame("events", inputfile) +d = d.Define("multiplicityS", "SimCaloHit_energy.size()")\ + .Define("multiplicityC", "SimCaloHitC_energy.size()")\ + .Define("etaHitsS", "SimCaloHit_eta.size()")\ + .Define("phiHitsS", "SimCaloHit_phi.size()")\ + .Define("enetotF", "Sum(SimCaloHit_energy[SimCaloHit_depth == 1])")\ + .Define("enetotR", "Sum(SimCaloHit_energy[SimCaloHit_depth == 2])")\ + .Define("enetot", "Sum(SimCaloHit_energy)")\ + .Define("HitEnergyMax", "Max(SimCaloHit_energy)")\ + .Define("indexHitEnergyMax", "ArgMax(SimCaloHit_energy)")\ + .Define("etaMaxHit", "SimCaloHit_eta[indexHitEnergyMax]")\ + .Define("phiMaxHit", "SimCaloHit_phi[indexHitEnergyMax]")\ + .Define("rMaxHit", "SimCaloHit_r[indexHitEnergyMax]")\ + .Define("xMaxHit", "SimCaloHit_x[indexHitEnergyMax]")\ + .Define("yMaxHit", "SimCaloHit_y[indexHitEnergyMax]")\ + .Define("indexHitCherMax", "ArgMax(SimCaloHitC_energy)")\ + .Define("etaMaxHitC", "SimCaloHit_eta[indexHitCherMax]")\ + .Define("phiMaxHitC", "SimCaloHit_phi[indexHitCherMax]")\ + .Define("rMaxHitC", "SimCaloHit_r[indexHitCherMax]")\ + .Define("xMaxHitC", "SimCaloHit_x[indexHitCherMax]")\ + .Define("yMaxHitC", "SimCaloHit_y[indexHitCherMax]") + + + +#float S_F = gRandom->Poisson(eneF*LO)/this_ene; +# float S_R = gRandom->Poisson(eneR*LO)/this_ene; +# float S = S_F+S_R; +c = ROOT.TCanvas("c","",800,600) + +#Get the En from the file +pGunEnergy = d.Histo1D ( ("pGunEnergy", "pGunEnergy", 50, 0, 50) , "MCparticle_energy" ).GetMean() + +#Define the histos 1D + + +myhistos1D = {} +#myhistos1D.append(d.Histo1D ( , "" )) + +myhistos1D["hMultiplicityS"] = d.Histo1D ( ("hMultiplicityS", "hMultiplicityS", 1500, 0, 1500) , "multiplicityS" ) +myhistos1D["hMultiplicityC"] = d.Histo1D ( ("hMultiplicityC", "hMultiplicityC", 1500, 0, 1500) , "multiplicityC" ) +myhistos1D["hEtaHits"] = d.Histo1D (("hEtaHits", "hEtaHits", 100, -10, 10) , "etaHitsS" ) +myhistos1D["hPhiHits"] = d.Histo1D (("hPhiHits", "hPhiHits", 100, -6.28, 6.28), "phiHitsS" ) +myhistos1D["hEneHits"] = d.Histo1D (("hEneHits", ";Hit energy [GeV]; Counts", 1000, 0, pGunEnergy*1.25), "SimCaloHit_energy" ) +myhistos1D["hTotEne"] = d.Histo1D (("hTotEne", " ;Total energy deposited [GeV]; Counts", 1000, 0, pGunEnergy*1.25), "enetot" ) +myhistos1D["hTotEneF"] = d.Histo1D (("hTotEneF"," ;Total energy deposited [GeV]; Counts", 1000, 0, pGunEnergy*1.25), "enetotF" ) +myhistos1D["hTotEneR"] = d.Histo1D (("hTotEneR"," ;Total energy deposited [GeV]; Counts", 1000, 0, pGunEnergy*1.25), "enetotR" ) + + +#Define the histos 2D +myhistos2D = {} +#myhistos2D.append(d.Histo2D (, "" )) +myhistos2D["hScatterMultiplicity"] = d.Histo2D (("hScatterMultiplicity", ";S hits multiplicity; C hits multiplicity", 300, 0, 1500, 300, 0, 1500), "multiplicityS", "multiplicityC" ) +#myhistos2D.append(d.Histo2D (("hScatterCS", "", 1000, 0, 20, 1000, 0, 1000), "" , "")) +myhistos2D["hTotEne_vs_eta"] = d.Histo2D (("hTotEne_vs_eta", ";#eta;Tot energy [GeV]", 100, -3.2, 3.2 , 1000, 0, pGunEnergy*1.25), "etaMaxHit" , "enetot") +myhistos2D["hTotEne_vs_phi"] = d.Histo2D (("hTotEne_vs_phi", ";#Phi;Tot energy [GeV]", 100, -3.2, 3.2 , 1000, 0, pGunEnergy*1.25), "phiMaxHit" , "enetot" ) + +myprofs = [] +#Define the profile +#myprofs.append(d.Profile1D(("hEne_vs_eta", "hEne_vs_eta", 100, -10, 10), "", "")) + +### Now drawing +ROOT.gStyle.SetOptStat(1) + +for key, h in myhistos1D.items(): + c.Clear() + c.cd() + + #h.GetXaxis().SetRangeUser(0, h.GetMean()+ h.GetRMS()*5) + h.Draw() + c.SetLogy(0) + c.SaveAs("%s/histo1D_%s.png"%(outdir,key)) + + c.SetLogy() + c.SaveAs("%s/log_histo1D_%s.png"%(outdir,key)) + +ROOT.gStyle.SetOptStat(0) + +for key, h in myhistos2D.items(): + c.Clear() + c.cd() + h.Draw("COLZ") + + c.SetLogy(0) + c.SaveAs("%s/histo2D_%s.png"%(outdir,key)) + + +### Ene Hits +myhistos1D["hEneHits"].SetStats(0) +c.SetLogy() +c.SetLogx() +myhistos1D["hEneHits"].GetXaxis().SetRangeUser(0.05,pGunEnergy*1.1) +myhistos1D["hEneHits"].Draw("") + +c.SaveAs("%s/cSCEP_EneHits.png"%outdir) + + +c.SetLogx(0) +### Energy +myhistos1D["hTotEne"].SetStats(0) + +ROOT.gStyle.SetOptStat(1) +myhistos1D["hTotEne"].SetStats(1) +myhistos1D["hTotEne"].GetXaxis().SetRangeUser(pGunEnergy*0.6,pGunEnergy*1.1) +myhistos1D["hTotEne"].Draw("") + +c.SetLogy(1) +c.SaveAs("%s/cSCEP_TotEnergy.png"%outdir) + + +c.Clear() +c.SetLogy(0) + +### Histo multiplicity +ROOT.gStyle.SetOptStat(0) +myhistos2D["hScatterMultiplicity"].SetStats(0) +myhistos2D["hScatterMultiplicity"].GetXaxis().SetRangeUser(10, myhistos1D["hMultiplicityS"].GetMean()+ myhistos1D["hMultiplicityS"].GetRMS()*5) +myhistos2D["hScatterMultiplicity"].GetYaxis().SetRangeUser(10, myhistos1D["hMultiplicityS"].GetMean()+ myhistos1D["hMultiplicityS"].GetRMS()*5) +#myhistos2D["hScatterMultiplicity"].GetYaxis().SetRangeUser(0, myhistos1D["hMultiplicityC"].GetMean()+ myhistos1D["hMultiplicityC"].GetRMS()*5) +myhistos2D["hScatterMultiplicity"].Draw("COLZ") +#c.SetLogz() +c.SaveAs("%s/cSCEP_ScatterMultiplicity.png"%outdir) + +# ETA +ROOT.gStyle.SetOptStat(0) +myhistos2D["hTotEne_vs_eta"].SetStats(0) +myhistos2D["hTotEne_vs_eta"].GetYaxis().SetRangeUser(pGunEnergy*0.6,pGunEnergy*1.1) +myhistos2D["hTotEne_vs_eta"].Draw("COLZ") +#c.SetLogz() +c.SaveAs("%s/cSCEP_hTotEne_vs_eta.png"%outdir) + +#PHI +myhistos2D["hTotEne_vs_phi"].SetStats(0) +myhistos2D["hTotEne_vs_phi"].GetYaxis().SetRangeUser(pGunEnergy*0.6,pGunEnergy*1.1) +myhistos2D["hTotEne_vs_phi"].Draw("COLZ") +#c.SetLogz() +c.SaveAs("%s/cSCEP_hTotEne_vs_phi.png"%outdir) + + + +##### Energy sharing plots +c.Clear() + +ROOT.gStyle.SetOptStat(0) +myhistos1D["hTotEne"].Draw() +myhistos1D["hTotEne"].GetXaxis().SetRangeUser(0,pGunEnergy*1.1) +myhistos1D["hTotEne"].SetStats(0) +#myhistos1D["hTotEne"].SetTitle(0) +myhistos1D["hTotEne"].SetLineWidth(2) +myhistos1D["hTotEne"].GetXaxis().SetTitle("Energy deposited in SCEPCAL [GeV]") +myhistos1D["hTotEne"].GetYaxis().SetTitle("Counts") +myhistos1D["hTotEne"].GetYaxis().SetRangeUser(1, myhistos1D["hTotEne"].GetMaximum()*5) +myhistos1D["hTotEne"].SetLineColor(1) + + + +myhistos1D["hTotEneF"].SetStats(0) +myhistos1D["hTotEneF"].SetLineColor(416 + 1) +myhistos1D["hTotEneF"].SetLineWidth(2) +myhistos1D["hTotEneF"].Draw("same") + + +myhistos1D["hTotEneR"].SetStats(0) +myhistos1D["hTotEneR"].SetLineColor(600) +myhistos1D["hTotEneR"].SetLineWidth(2) +myhistos1D["hTotEneR"].Draw("same") + +leg = ROOT.TLegend(0.15,0.68,0.45,0.88) + +h = myhistos1D["hTotEne"] +leg.AddEntry( "hTotEneF" , "Front crystal", "lp") +leg.AddEntry("hTotEneR" , "Rear crystal", "lp") +leg.AddEntry( "hTotEne" , "Total", "lp") + +leg.Draw("same") + +c.SetLogy() +c.SaveAs("%s/cSCEP_EneSharing.png"%outdir) + +### nHits plots +c.Clear() + +ROOT.gStyle.SetOptStat(0) +myhistos1D["hMultiplicityS"].Draw("") +myhistos1D["hMultiplicityS"].SetStats(0) +myhistos1D["hMultiplicityC"].SetStats(0) +myhistos1D["hMultiplicityS"].GetXaxis().SetTitle("Hits Multiplicity") +myhistos1D["hMultiplicityS"].GetYaxis().SetTitle("Counts") +myhistos1D["hMultiplicityS"].GetYaxis().SetRangeUser(1, myhistos1D["hTotEne"].GetMaximum()*5) + +myhistos1D["hMultiplicityS"].SetLineColor(416 + 1) +myhistos1D["hMultiplicityS"].SetLineWidth(2) +myhistos1D["hMultiplicityS"].Draw("") + +myhistos1D["hMultiplicityC"].SetLineColor(600) +myhistos1D["hMultiplicityC"].SetLineWidth(2) +myhistos1D["hMultiplicityC"].Draw("same") + + +leg = ROOT.TLegend(0.15,0.68,0.45,0.88) + +leg.AddEntry( "hMultiplicityS" , "Scintillation", "lp") +leg.AddEntry( "hMultiplicityC" , "Cherenkov", "lp") + +leg.Draw("same") + +c.SetLogy() +c.SaveAs("%s/cSCEP_Multiplicity.png"%outdir) + + + + + + + From fb93ff7cc6523bcd29b628a9c1c5b0033260ccc2 Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Thu, 16 Nov 2023 18:07:59 +0100 Subject: [PATCH 5/7] minor fix --- SCEPCal_plots/ntuplizer.py | 3 -- SCEPCal_plots/simplePlotter.py | 14 ++------- .../dataframe/FCCAnalyses/CaloNtupleizer.h | 4 --- analyzers/dataframe/src/CaloNtupleizer.cc | 30 +------------------ 4 files changed, 3 insertions(+), 48 deletions(-) diff --git a/SCEPCal_plots/ntuplizer.py b/SCEPCal_plots/ntuplizer.py index f280f9a0e8..3a6f0d2b9e 100644 --- a/SCEPCal_plots/ntuplizer.py +++ b/SCEPCal_plots/ntuplizer.py @@ -10,9 +10,6 @@ _fcc = ROOT.dummyLoader ROOT.gInterpreter.Declare("using namespace FCCAnalyses;") -#to add layer info - not working -#ROOT.CaloNtupleizer.loadGeometry('/afs/cern.ch/work/f/fcetorel/private/work2/SCEPCAL/Detector/SCEPCAL/compact/SCEPCAL.xml', 'SCEPCALreadout') - #/afs/cern.ch/work/f/fcetorel/private/work2/SCEPCAL/SCEPCALsim/SCEPCALsimG4Components/test/data/ parser = argparse.ArgumentParser() diff --git a/SCEPCal_plots/simplePlotter.py b/SCEPCal_plots/simplePlotter.py index 02168a2723..f269aa092d 100755 --- a/SCEPCal_plots/simplePlotter.py +++ b/SCEPCal_plots/simplePlotter.py @@ -26,12 +26,10 @@ args = parser.parse_args() # input file should be flat ntupla from ntuplizer -#outdir = "/eos/user/f/fcetorel/www/SCEPCal_Sim/test_plots/newCommitsChecks/" -#inputfile = "flatNtupla_scepcal_100evs_newPull.root" outdir = args.outFolder inputfile = args.inputFile os.makedirs(outdir, exist_ok=True) -os.system("cp /eos/user/f/fcetorel/www/index.php %s"%outdir) +#os.system("cp /eos/user/f/fcetorel/www/index.php %s"%outdir) @@ -60,9 +58,6 @@ -#float S_F = gRandom->Poisson(eneF*LO)/this_ene; -# float S_R = gRandom->Poisson(eneR*LO)/this_ene; -# float S = S_F+S_R; c = ROOT.TCanvas("c","",800,600) #Get the En from the file @@ -92,10 +87,6 @@ myhistos2D["hTotEne_vs_eta"] = d.Histo2D (("hTotEne_vs_eta", ";#eta;Tot energy [GeV]", 100, -3.2, 3.2 , 1000, 0, pGunEnergy*1.25), "etaMaxHit" , "enetot") myhistos2D["hTotEne_vs_phi"] = d.Histo2D (("hTotEne_vs_phi", ";#Phi;Tot energy [GeV]", 100, -3.2, 3.2 , 1000, 0, pGunEnergy*1.25), "phiMaxHit" , "enetot" ) -myprofs = [] -#Define the profile -#myprofs.append(d.Profile1D(("hEne_vs_eta", "hEne_vs_eta", 100, -10, 10), "", "")) - ### Now drawing ROOT.gStyle.SetOptStat(1) @@ -138,7 +129,7 @@ ROOT.gStyle.SetOptStat(1) myhistos1D["hTotEne"].SetStats(1) -myhistos1D["hTotEne"].GetXaxis().SetRangeUser(pGunEnergy*0.6,pGunEnergy*1.1) +myhistos1D["hTotEne"].GetXaxis().SetRangeUser(pGunEnergy*0.7,pGunEnergy*1.1) myhistos1D["hTotEne"].Draw("") c.SetLogy(1) @@ -153,7 +144,6 @@ myhistos2D["hScatterMultiplicity"].SetStats(0) myhistos2D["hScatterMultiplicity"].GetXaxis().SetRangeUser(10, myhistos1D["hMultiplicityS"].GetMean()+ myhistos1D["hMultiplicityS"].GetRMS()*5) myhistos2D["hScatterMultiplicity"].GetYaxis().SetRangeUser(10, myhistos1D["hMultiplicityS"].GetMean()+ myhistos1D["hMultiplicityS"].GetRMS()*5) -#myhistos2D["hScatterMultiplicity"].GetYaxis().SetRangeUser(0, myhistos1D["hMultiplicityC"].GetMean()+ myhistos1D["hMultiplicityC"].GetRMS()*5) myhistos2D["hScatterMultiplicity"].Draw("COLZ") #c.SetLogz() c.SaveAs("%s/cSCEP_ScatterMultiplicity.png"%outdir) diff --git a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h index 909b7ac70e..0a41b5d632 100644 --- a/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h +++ b/analyzers/dataframe/FCCAnalyses/CaloNtupleizer.h @@ -39,17 +39,13 @@ ROOT::VecOps::RVec getSimCaloHit_x (const ROOT::VecOps::RVec getSimCaloHit_y (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_z (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_phi (const ROOT::VecOps::RVec& in); -ROOT::VecOps::RVec getSimCaloHit_phiBin (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_theta (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec& in); -ROOT::VecOps::RVec getSimCaloHit_etaBin (const ROOT::VecOps::RVec& in); -ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_depth (const ROOT::VecOps::RVec& in,const int decodingVal); ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps::RVec& in); - // calo hits (single cells) ROOT::VecOps::RVec getCaloHit_x (const ROOT::VecOps::RVec& in); ROOT::VecOps::RVec getCaloHit_y (const ROOT::VecOps::RVec& in); diff --git a/analyzers/dataframe/src/CaloNtupleizer.cc b/analyzers/dataframe/src/CaloNtupleizer.cc index 14d5ad0769..8996c08e64 100644 --- a/analyzers/dataframe/src/CaloNtupleizer.cc +++ b/analyzers/dataframe/src/CaloNtupleizer.cc @@ -76,15 +76,6 @@ ROOT::VecOps::RVec getSimCaloHit_phi (const ROOT::VecOps::RVec getSimCaloHit_phiBin (const ROOT::VecOps::RVec& in){ -// ROOT::VecOps::RVec result; -// for (auto & p: in){ -// dd4hep::DDSegmentation::CellID cellId = p.cellID; -// result.push_back(m_decoder->getSim(cellId, "phi")); -// } -// return result; -//} - ROOT::VecOps::RVec getSimCaloHit_theta (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -105,15 +96,6 @@ ROOT::VecOps::RVec getSimCaloHit_eta (const ROOT::VecOps::RVec getSimCaloHit_etaBin (const ROOT::VecOps::RVec& in){ -// ROOT::VecOps::RVec result; -// for (auto & p: in){ -// dd4hep::DDSegmentation::CellID cellId = p.cellID; -// result.push_back(m_decoder->get(cellId, "eta")); -// } -// return result; -//} - ROOT::VecOps::RVec getSimCellID (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -122,8 +104,6 @@ ROOT::VecOps::RVec getSimCellID (const ROOT::VecOps::RVec getSimCaloHit_energy (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -140,15 +120,6 @@ ROOT::VecOps::RVec getSimCaloHit_depth (const ROOT::VecOps::RVec getSimCaloHit_layer (const ROOT::VecOps::RVec& in){ - ROOT::VecOps::RVec result; - for (auto & p: in){ - dd4hep::DDSegmentation::CellID cellId = p.cellID; - result.push_back(m_decoder->get(cellId, "layer")); - } - return result; -} - ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; for (auto & p: in){ @@ -158,6 +129,7 @@ ROOT::VecOps::RVec getSimCaloHit_positionVector3 (const ROOT::VecOps:: } return result; } + // calo hit ROOT::VecOps::RVec getCaloHit_x (const ROOT::VecOps::RVec& in){ ROOT::VecOps::RVec result; From e2216c65ed55840041d4230673d24036dc09455d Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Wed, 6 Dec 2023 16:24:27 +0100 Subject: [PATCH 6/7] moved SCEPCal_plots to examples/FCCee --- {SCEPCal_plots => examples/FCCee/SCEPCal_plots}/ntuplizer.py | 0 {SCEPCal_plots => examples/FCCee/SCEPCal_plots}/simplePlotter.py | 0 2 files changed, 0 insertions(+), 0 deletions(-) rename {SCEPCal_plots => examples/FCCee/SCEPCal_plots}/ntuplizer.py (100%) rename {SCEPCal_plots => examples/FCCee/SCEPCal_plots}/simplePlotter.py (100%) diff --git a/SCEPCal_plots/ntuplizer.py b/examples/FCCee/SCEPCal_plots/ntuplizer.py similarity index 100% rename from SCEPCal_plots/ntuplizer.py rename to examples/FCCee/SCEPCal_plots/ntuplizer.py diff --git a/SCEPCal_plots/simplePlotter.py b/examples/FCCee/SCEPCal_plots/simplePlotter.py similarity index 100% rename from SCEPCal_plots/simplePlotter.py rename to examples/FCCee/SCEPCal_plots/simplePlotter.py From a693dd53d4e1a2a7006d15781318b2c8299b664b Mon Sep 17 00:00:00 2001 From: Flavia Cetorelli Date: Wed, 6 Dec 2023 16:32:33 +0100 Subject: [PATCH 7/7] added option to specify output dir --- examples/FCCee/SCEPCal_plots/ntuplizer.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/FCCee/SCEPCal_plots/ntuplizer.py b/examples/FCCee/SCEPCal_plots/ntuplizer.py index 3a6f0d2b9e..8e601fc50a 100644 --- a/examples/FCCee/SCEPCal_plots/ntuplizer.py +++ b/examples/FCCee/SCEPCal_plots/ntuplizer.py @@ -10,11 +10,11 @@ _fcc = ROOT.dummyLoader ROOT.gInterpreter.Declare("using namespace FCCAnalyses;") -#/afs/cern.ch/work/f/fcetorel/private/work2/SCEPCAL/SCEPCALsim/SCEPCALsimG4Components/test/data/ parser = argparse.ArgumentParser() parser.add_argument("-inputFilesName", help = "name of the input rootfiles", type = str) -parser.add_argument("-baseFolder",help = "folder for the rootfiles", type = str) +parser.add_argument("-baseFolder",help = "input folder", type = str) +parser.add_argument("-outFolder",help = "output folder", type = str) args = parser.parse_args() NBITS = 27 @@ -50,7 +50,7 @@ ) -outfilename = 'flatNtupla_'+ args.inputFilesName +outfilename = args.outFolder+'/flatNtupla_'+ args.inputFilesName print(f'Writing snapshot to disk ... \t{outfilename}') df.Snapshot('events', outfilename,