From 023b1933e3d3e8dea7a30d626db2f9cc3b59d1fc Mon Sep 17 00:00:00 2001 From: Edi Date: Thu, 10 Mar 2022 17:00:06 +0100 Subject: [PATCH 1/6] updated ghostMatching defn --- analyzers/dataframe/JetTaggingUtils.cc | 204 +++++++++++++++++++++++++ analyzers/dataframe/JetTaggingUtils.h | 38 +++++ 2 files changed, 242 insertions(+) diff --git a/analyzers/dataframe/JetTaggingUtils.cc b/analyzers/dataframe/JetTaggingUtils.cc index 37c497e518..c46c36f643 100644 --- a/analyzers/dataframe/JetTaggingUtils.cc +++ b/analyzers/dataframe/JetTaggingUtils.cc @@ -61,6 +61,210 @@ JetTaggingUtils::get_flavour(ROOT::VecOps::RVec in, return result; } + +get_ghostFlavour::get_ghostFlavour(int algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, float arg_add1, float arg_add2) +{m_algo = algo; m_radius = arg_radius; m_exclusive = arg_exclusive; m_cut = arg_cut; m_sorted = arg_sorted; m_recombination = arg_recombination; m_add1 = arg_add1; m_add2 = arg_add2;} + +ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec Particle, ROOT::VecOps::RVec ind, std::vector pseudoJets, int partonFlag) { + + + + ghostFlavour result; + + unsigned index = pseudoJets.size(); + std::vector pdg(pseudoJets.size(),0); + std::vector ghostStatus(pseudoJets.size(),0); + std::vector MCindex(pseudoJets.size(),-1); + + + std::vector partonStatus = {20, 30}; + + // In below loop ghosts are selected from MCParticle collection + for (size_t i = 0; i < Particle.size(); ++i) { + bool isGhost = false; + + if(!partonFlag){ + // Ghost partons as any partons that do not have partons as daughters + if (std::abs(Particle[i].PDG)<=5||Particle[i].PDG==21) { + isGhost = true; + auto daughters = MCParticle::get_list_of_particles_from_decay(i, Particle, ind); + for(auto& daughter_index : daughters){ + if (std::abs(Particle[daughter_index].PDG)<=5||Particle[daughter_index].PDG==21) { + isGhost = false; + break; + } + } + if (isGhost) ghostStatus.push_back(1); + } + } + else{ + // Ghost partons are selected via chosen Pythia8 status codes + if ((Particle[i].generatorStatuspartonStatus[0])) { + isGhost = true; + ghostStatus.push_back(1); + } + } + + + // Ghost hadrons are selected as b/c hadrons that do not have b/c hadrons as daughters + if ((std::abs(int((Particle[i].PDG/100))%10)==5)||(std::abs(int((Particle[i].PDG/1000))%10)==5)){ + isGhost = true; + auto daughters = MCParticle::get_list_of_particles_from_decay(i, Particle, ind); + for(auto& daughter_index : daughters){ + if((std::abs(int((Particle[daughter_index].PDG/100))%10)==5)||(std::abs(int((Particle[daughter_index].PDG/1000))%10)==5)){ + isGhost = false; + break; + } + } + if (isGhost) ghostStatus.push_back(2); + + } + else if ((std::abs(int((Particle[i].PDG/100))%10)==4)||(std::abs(int((Particle[i].PDG/1000))%10)==4)){ + isGhost = true; + auto daughters = MCParticle::get_list_of_particles_from_decay(i, Particle, ind); + for(auto& daughter_index : daughters){ + if((std::abs(int((Particle[daughter_index].PDG/100))%10)==4)||(std::abs(int((Particle[daughter_index].PDG/1000))%10)==4)){ + isGhost = false; + break; + } + } + if (isGhost) ghostStatus.push_back(2); + } + + // Ghosts 4-mom is scaled by 10^-18 + if (isGhost){ + pdg.push_back(Particle[i].PDG); + MCindex.push_back(i); + // the double conversion here is verbose but for precision if I recall... + double px = Particle[i].momentum.x;//*pow(10, -1); + double py = Particle[i].momentum.y; + double pz = Particle[i].momentum.z; + double m = Particle[i].mass; + double E = sqrt(px*px + py*py + pz*pz + m*m); + pseudoJets.emplace_back(px*pow(10, -18), py*pow(10, -18), pz*pow(10, -18), E*pow(10, -18)); + pseudoJets.back().set_user_index(index); + ++index; + } + } + + result.ghostStatus = ghostStatus; + result.MCindex = MCindex; + + + // Jet clustering algorithm is run according to user choice m_algo + JetClusteringUtils::FCCAnalysesJet FCCAnalysesGhostJets; + + if (m_algo == 0) FCCAnalysesGhostJets = JetClustering::clustering_kt(m_radius, m_exclusive, m_cut, m_sorted, m_recombination)(pseudoJets); + else if (m_algo == 1) FCCAnalysesGhostJets = JetClustering::clustering_antikt(m_radius, m_exclusive, m_cut, m_sorted, m_recombination)(pseudoJets); + else if (m_algo == 2) FCCAnalysesGhostJets = JetClustering::clustering_cambridge(m_radius, m_exclusive, m_cut, m_sorted, m_recombination)(pseudoJets); + else if (m_algo == 3) FCCAnalysesGhostJets = JetClustering::clustering_ee_kt(m_exclusive, m_cut, m_sorted, m_recombination)(pseudoJets); + else if (m_algo == 4) FCCAnalysesGhostJets = JetClustering::clustering_ee_genkt(m_radius, m_exclusive, m_cut, m_sorted, m_recombination, m_add1)(pseudoJets); + else if (m_algo == 5) FCCAnalysesGhostJets = JetClustering::clustering_genkt(m_radius, m_exclusive, m_cut, m_sorted, m_recombination, m_add1)(pseudoJets); + else if (m_algo == 6) FCCAnalysesGhostJets = JetClustering::clustering_valencia(m_radius, m_exclusive, m_cut, m_sorted, m_recombination, m_add1, m_add2)(pseudoJets); + else if (m_algo == 7) FCCAnalysesGhostJets = JetClustering::clustering_jade(m_radius, m_exclusive, m_cut, m_sorted, m_recombination)(pseudoJets); + else return result; + + + result.jets = FCCAnalysesGhostJets; + + // Jet constituents and pseudojets are read from resulting jet struct + auto ghostJets_ee_kt = JetClusteringUtils::get_pseudoJets(FCCAnalysesGhostJets); + + auto jetconstituents = JetClusteringUtils::get_constituents(FCCAnalysesGhostJets); + + + + + + + // Flav vector is defined before jets are check for clustered ghosts + std::vector> flav_vector; + + + std::vector partonFlavs; + std::vector hadronFlavs; + for (auto& consti_index : jetconstituents) { + + float partonFlav = 0; + float partonMom2 = 0; + float partonMom2_b = 0; + float partonMom2_c = 0; + float hadronFlav = 0; + float hadronMom2_b = 0; + float hadronMom2_c = 0; + + for (auto& i : consti_index) { + + //Parton-flav loop + if (ghostStatus[i]==1) { + if (std::abs(pdg[i])==5){ + if (pseudoJets[i].modp2()>partonMom2_b){ + partonFlav = pdg[i]; + partonMom2_b = pseudoJets[i].modp2(); + } + } + else if ((std::abs(pdg[i])==4) && (std::abs(partonFlav)<5)) { + if (pseudoJets[i].modp2()>partonMom2_c){ + partonFlav = pdg[i]; + partonMom2_c = pseudoJets[i].modp2(); + } + } + else if (((std::abs(pdg[i])==3) || (std::abs(pdg[i])==2) || (std::abs(pdg[i])==1) || (std::abs(pdg[i])==21)) && ((std::abs(partonFlav)<4) || (std::abs(partonFlav)==21))) { + if (pseudoJets[i].modp2()>partonMom2){ + partonFlav = pdg[i]; + partonMom2 = pseudoJets[i].modp2(); + } + } + } + + // Hadron-flav loop + if (ghostStatus[i]==2) { + if ((std::abs(int((pdg[i]/100))%10)==5)||(std::abs(int((pdg[i]/1000))%10)==5)){ + if (pseudoJets[i].modp2()>hadronMom2_b){ + hadronFlav = ((pdg[i]<0)-(pdg[i]>0))*5; + hadronMom2_b = pseudoJets[i].modp2(); + } + } + else if (((std::abs(int((pdg[i]/100))%10)==4)||(std::abs(int((pdg[i]/1000))%10)==4)) && (std::abs(hadronFlav)<5)) { + if (pseudoJets[i].modp2()>hadronMom2_c){ + hadronFlav = ((pdg[i]>0)-(pdg[i]<0))*4; + hadronMom2_c = pseudoJets[i].modp2(); + } + } + } + } + partonFlavs.push_back(partonFlav); + hadronFlavs.push_back(hadronFlav); + } + + + flav_vector.push_back(partonFlavs); + flav_vector.push_back(hadronFlavs); + + result.flavour = flav_vector; + return result; + +} + +std::vector> JetTaggingUtils::get_flavour(ghostFlavour ghostStruct){ + return ghostStruct.flavour; +} + +JetClusteringUtils::FCCAnalysesJet JetTaggingUtils::get_jets(ghostFlavour ghostStruct){ + return ghostStruct.jets; +} + +std::vector JetTaggingUtils::get_ghostStatus(ghostFlavour ghostStruct){ + return ghostStruct.ghostStatus; +} + +std::vector JetTaggingUtils::get_MCindex(ghostFlavour ghostStruct){ + return ghostStruct.MCindex; +} + + + + ROOT::VecOps::RVec JetTaggingUtils::get_btag(ROOT::VecOps::RVec in, float efficiency, float mistag_c, diff --git a/analyzers/dataframe/JetTaggingUtils.h b/analyzers/dataframe/JetTaggingUtils.h index 30bf4b8d58..09809c8a5c 100644 --- a/analyzers/dataframe/JetTaggingUtils.h +++ b/analyzers/dataframe/JetTaggingUtils.h @@ -20,6 +20,44 @@ namespace JetTaggingUtils{ //Get flavour association of jet ROOT::VecOps::RVec get_flavour(ROOT::VecOps::RVec in, ROOT::VecOps::RVec MCin); + + struct ghostFlavour { + std::vector> flavour; + JetClusteringUtils::FCCAnalysesJet jets; + std::vector ghostStatus; + std::vector MCindex; + }; + + + //Get ghost flavour (MC flavour) of jet described here: .. + //a struct is returned containing + // - vector of flavours: <, > + // - FCCAnalysesJet struct containing resulting jets (on which functions defined in JetClusteringUtils.cc can be called) + // - vector of jet constituent indices <, , ...> + // - ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) + // - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) + struct get_ghostFlavour { + get_ghostFlavour(int algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, float arg_add1=0, float arg_add2=0); + + int m_algo = 0; ///< flag to select jet clustering algorithm defined in JetClustering.cc (0 = kt, 1 = antikt, 2 = cambridge, 3 = eekt, 4 = ee genkt, 5 = genkt, 6 = valencia, 7 = jade) + float m_radius = 0.5; ///< jet cone radius (note that this variable should be passed even when not used e.g. for eekt) + int m_exclusive = 0; ///< flag for exclusive jet clustering. Possible choices are 0=inclusive clustering, 1=exclusive clustering that would be obtained when running the algorithm with the given dcut, 2=exclusive clustering when the event is clustered (in the exclusive sense) to exactly njets, 3=exclusive clustering when the event is clustered (in the exclusive sense) up to exactly njets, 4=exclusive jets obtained at the given ycut + float m_cut = 5.; ///< pT cut for m_exclusive=0, dcut for m_exclusive=1, N jets for m_exlusive=2, N jets for m_exclusive=3, ycut for m_exclusive=4 + int m_sorted = 0; ///< pT ordering=0, E ordering=1 + int m_recombination = 0; ///< E_scheme=0, pt_scheme=1, pt2_scheme=2, Et_scheme=3, Et2_scheme=4, BIpt_scheme=5, BIpt2_scheme=6, E0_scheme=10, p_scheme=11 + float m_add1 = 0.; /// first additional parameter (to be used in selected clustering scheme) + float m_add2 = 0.; /// second additional parameter + ghostFlavour operator() (ROOT::VecOps::RVec Particle, ROOT::VecOps::RVec ind, std::vector pseudoJets, int partonFlag); + }; + + std::vector> get_flavour(ghostFlavour ghostStruct); + + JetClusteringUtils::FCCAnalysesJet get_jets(ghostFlavour ghostStruct); + + std::vector get_ghostStatus(ghostFlavour ghostStruct); + + std::vector get_MCindex(ghostFlavour ghostStruct); + //Get b-tags with an efficiency applied ROOT::VecOps::RVec get_btag(ROOT::VecOps::RVec in, float efficiency, float mistag_c=0., float mistag_l=0., float mistag_g=0.); //Get c-tags with an efficiency applied From 23a6cf2dab6f149536e6998017450b8c1cd5a89e Mon Sep 17 00:00:00 2001 From: Kunal Gautam Date: Thu, 10 Mar 2022 17:53:30 +0100 Subject: [PATCH 2/6] header update --- analyzers/dataframe/JetTaggingUtils.cc | 6 +++--- analyzers/dataframe/JetTaggingUtils.h | 4 ++++ 2 files changed, 7 insertions(+), 3 deletions(-) diff --git a/analyzers/dataframe/JetTaggingUtils.cc b/analyzers/dataframe/JetTaggingUtils.cc index c46c36f643..015a71b822 100644 --- a/analyzers/dataframe/JetTaggingUtils.cc +++ b/analyzers/dataframe/JetTaggingUtils.cc @@ -71,7 +71,7 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec pdg(pseudoJets.size(),0); std::vector ghostStatus(pseudoJets.size(),0); std::vector MCindex(pseudoJets.size(),-1); @@ -79,7 +79,7 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec partonStatus = {20, 30}; - // In below loop ghosts are selected from MCParticle collection + // In loop below ghosts are selected from MCParticle collection for (size_t i = 0; i < Particle.size(); ++i) { bool isGhost = false; @@ -177,7 +177,7 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec> flav_vector; diff --git a/analyzers/dataframe/JetTaggingUtils.h b/analyzers/dataframe/JetTaggingUtils.h index 09809c8a5c..20523bfbff 100644 --- a/analyzers/dataframe/JetTaggingUtils.h +++ b/analyzers/dataframe/JetTaggingUtils.h @@ -5,8 +5,12 @@ #include "Math/Vector4D.h" #include "ROOT/RVec.hxx" #include "edm4hep/MCParticleData.h" +#include "edm4hep/ReconstructedParticleData.h" #include "fastjet/JetDefinition.hh" #include "TRandom3.h" +#include "JetClustering.h" +#include "JetClusteringUtils.h" +#include "MCParticle.h" /** Jet tagging utilities interface. This represents a set functions and utilities to perfom jet tagging from a list of jets. From 6977d1f25d768e067b3647d8e4b63cab90f00e7f Mon Sep 17 00:00:00 2001 From: Eduardo Ploerer Date: Thu, 10 Mar 2022 18:37:24 +0100 Subject: [PATCH 3/6] ghost Flavour defn --- examples/FCCee/ghostFlavour/analysis.py | 170 ++++++++++++++++++++++++ examples/FCCee/ghostFlavour/preSel.py | 17 +++ 2 files changed, 187 insertions(+) create mode 100644 examples/FCCee/ghostFlavour/analysis.py create mode 100644 examples/FCCee/ghostFlavour/preSel.py diff --git a/examples/FCCee/ghostFlavour/analysis.py b/examples/FCCee/ghostFlavour/analysis.py new file mode 100644 index 0000000000..37becde23f --- /dev/null +++ b/examples/FCCee/ghostFlavour/analysis.py @@ -0,0 +1,170 @@ +import sys +import ROOT + +print ("Load cxx analyzers ... ",) +ROOT.gSystem.Load("libedm4hep") +ROOT.gSystem.Load("libpodio") +ROOT.gSystem.Load("libFCCAnalyses") +ROOT.gErrorIgnoreLevel = ROOT.kFatal +_edm = ROOT.edm4hep.ReconstructedParticleData() +_pod = ROOT.podio.ObjectID() +_fcc = ROOT.dummyLoader + +print ('edm4hep ',_edm) +print ('podio ',_pod) +print ('fccana ',_fcc) + +class analysis(): + + #__________________________________________________________ + def __init__(self, inputlist, outname, ncpu): + self.outname = outname + if ".root" not in outname: + self.outname+=".root" + + ROOT.ROOT.EnableImplicitMT(ncpu) + + self.df = ROOT.RDataFrame("events", inputlist) + print (" done") + #__________________________________________________________ + def run(self): + df2 = (self.df + #df2 = (self.df.Range(0, 10) + #alias for dealing with # in python + .Alias("Particle0", "Particle#0.index") + .Alias("Particle1", "Particle#1.index") + .Alias("Jet3","Jet#3.index") + + .Define("MC_px", "MCParticle::get_px(Particle)") + .Define("MC_py", "MCParticle::get_py(Particle)") + .Define("MC_pz", "MCParticle::get_pz(Particle)") + .Define("MC_p", "MCParticle::get_p(Particle)") + .Define("MC_e", "MCParticle::get_e(Particle)") + .Define("MC_m", "MCParticle::get_mass(Particle)") + .Define("MC_theta", "MCParticle::get_theta(Particle)") + .Define("MC_pdg", "MCParticle::get_pdg(Particle)") + .Define("MC_status", "MCParticle::get_genStatus(Particle)") + + + #define the RP px, py, pz and e + .Define("RP_px", "ReconstructedParticle::get_px(ReconstructedParticles)") + .Define("RP_py", "ReconstructedParticle::get_py(ReconstructedParticles)") + .Define("RP_pz", "ReconstructedParticle::get_pz(ReconstructedParticles)") + .Define("RP_m", "ReconstructedParticle::get_mass(ReconstructedParticles)") + + + + + ################ + #Jet Clustering# + ################ + + #build pseudo jets with the RP, using the interface that takes px,py,pz,m for better + #handling of rounding errors + .Define("pseudo_jets", "JetClusteringUtils::set_pseudoJets_xyzm(RP_px, RP_py, RP_pz, RP_m)") + + #EE-KT ALGORITHM + #run jet clustering with all MC particles. ee_kt_algorithm, exclusive clustering, exactly 2 jets, E-scheme + .Define("FCCAnalysesJets_ee_kt", "JetClustering::clustering_ee_kt(2, 2, 1, 0)(pseudo_jets)") + + #get the jets out of the structure + .Define("jets_ee_kt", "JetClusteringUtils::get_pseudoJets(FCCAnalysesJets_ee_kt)") + + #get the jet constituents out of the structure + .Define("jetconstituents_ee_kt", "JetClusteringUtils::get_constituents(FCCAnalysesJets_ee_kt)") + + #get some jet variables + .Define("jets_ee_kt_e", "JetClusteringUtils::get_e(jets_ee_kt)") + .Define("jets_ee_kt_px", "JetClusteringUtils::get_px(jets_ee_kt)") + .Define("jets_ee_kt_py", "JetClusteringUtils::get_py(jets_ee_kt)") + .Define("jets_ee_kt_pz", "JetClusteringUtils::get_pz(jets_ee_kt)") + .Define("jets_ee_kt_flavour", "JetTaggingUtils::get_flavour(jets_ee_kt, Particle)") + + ############### + #Ghost Flavour# + ############### + + #pseudo jets defined above, along with Particle, Particle 1 collections are passed to get_ghostflavour + .Define("ghostFlavour", "JetTaggingUtils::get_ghostFlavour(3, 0, 2, 2, 1, 0)(Particle, Particle1, pseudo_jets, 0)") + + #the flavour vector can be obtained from the struct by using get_flavour + .Define("ghostFlavour_flav", "JetTaggingUtils::get_flavour(ghostFlavour)") + + #an FCCAnalysesJet struct can be obtain by using get_jets + .Define("ghostFlavour_jets", "JetTaggingUtils::get_jets(ghostFlavour)") + #the pseudojets can then be obtained (as above) using get_pseudojets + .Define("ghostJets", "JetClusteringUtils::get_pseudoJets(ghostFlavour_jets)") + #with the kinematics of the pseudojets accessible via their own functions + .Define("ghostJets_e", "JetClusteringUtils::get_e(ghostJets)") + .Define("ghostJets_px", "JetClusteringUtils::get_px(ghostJets)") + .Define("ghostJets_py", "JetClusteringUtils::get_py(ghostJets)") + .Define("ghostJets_pz", "JetClusteringUtils::get_pz(ghostJets)") + + #the ghostStatus and MCindex vector can likewise be obtained + .Define("ghostFlavour_MCindex", "JetTaggingUtils::get_MCindex(ghostFlavour)") + .Define("ghostFlavour_ghostStatus", "JetTaggingUtils::get_ghostStatus(ghostFlavour)") + + + + + ) + + + + + # select branches for output file + branchList = ROOT.vector('string')() + for branchName in [ + "MC_px", + "MC_py", + "MC_pz", + "MC_p", + "MC_e", + "MC_theta", + "MC_pdg", + "MC_status", + + "jets_ee_kt_e", + "jets_ee_kt_px", + "jets_ee_kt_py", + "jets_ee_kt_pz", + "jets_ee_kt_flavour", + "jetconstituents_ee_kt", + + "ghostFlavour_flav", + "ghostFlavour_MCindex", + "ghostFlavour_ghostStatus", + + "ghostJets_e", + "ghostJets_px", + "ghostJets_py", + "ghostJets_pz", + + + ]: + branchList.push_back(branchName) + df2.Snapshot("events", self.outname, branchList) + +# example call for standalone file +# python examples/FCCee/top/hadronic/analysis.py /eos/experiment/fcc/ee/generation/DelphesEvents/fcc_tmp/p8_ee_tt_fullhad_ecm365/events_196309147.root +if __name__ == "__main__": + + if len(sys.argv)==1: + print ("usage:") + print ("python ",sys.argv[0]," file.root") + sys.exit(3) + infile = sys.argv[1] + outDir = sys.argv[0].replace(sys.argv[0].split('/')[0],'outputs/FCCee').replace('analysis.py','')+'/' + import os + os.system("mkdir -p {}".format(outDir)) + outfile = outDir+infile.split('/')[-1] + ncpus = 2 + analysis = analysis(infile, outfile, ncpus) + analysis.run() + + tf = ROOT.TFile(infile) + entries = tf.events.GetEntries() + p = ROOT.TParameter(int)( "eventsProcessed", entries) + outf=ROOT.TFile(outfile,"UPDATE") + p.Write() + diff --git a/examples/FCCee/ghostFlavour/preSel.py b/examples/FCCee/ghostFlavour/preSel.py new file mode 100644 index 0000000000..4bf1f32a89 --- /dev/null +++ b/examples/FCCee/ghostFlavour/preSel.py @@ -0,0 +1,17 @@ + +from config.common_defaults import deffccdicts +import os + +basedir=os.path.join(os.getenv('FCCDICTSDIR', deffccdicts), '') + "yaml/FCCee/spring2021/IDEA/" +print('basedir: '+str(basedir)) +outdir="outputs/FCCee/ghostFlavour" + +import multiprocessing +NUM_CPUS = int(multiprocessing.cpu_count()-2) + +process_list=['p8_ee_Zbb_ecm91'] +fraction=0.0000000001 + +import config.runDataFrame as rdf +myana=rdf.runDataFrame(basedir,process_list) +myana.run(ncpu=NUM_CPUS,fraction=fraction,outDir=outdir) From e15759df90217f05a9f9d85d79ea09509529cac6 Mon Sep 17 00:00:00 2001 From: Kunal Gautam Date: Mon, 14 Mar 2022 17:48:42 +0100 Subject: [PATCH 4/6] suggested changes + added jet constituents --- analyzers/dataframe/JetTaggingUtils.cc | 18 ++++++----- analyzers/dataframe/JetTaggingUtils.h | 42 ++++++++++++++++++-------- 2 files changed, 40 insertions(+), 20 deletions(-) diff --git a/analyzers/dataframe/JetTaggingUtils.cc b/analyzers/dataframe/JetTaggingUtils.cc index 015a71b822..031ec210b0 100644 --- a/analyzers/dataframe/JetTaggingUtils.cc +++ b/analyzers/dataframe/JetTaggingUtils.cc @@ -72,12 +72,12 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec pdg(pseudoJets.size(),0); - std::vector ghostStatus(pseudoJets.size(),0); - std::vector MCindex(pseudoJets.size(),-1); + ROOT::VecOps::RVec pdg(pseudoJets.size(),0); + ROOT::VecOps::RVec ghostStatus(pseudoJets.size(),0); + ROOT::VecOps::RVec MCindex(pseudoJets.size(),-1); - std::vector partonStatus = {20, 30}; + ROOT::VecOps::RVec partonStatus = {20, 30}; // In loop below ghosts are selected from MCParticle collection for (size_t i = 0; i < Particle.size(); ++i) { @@ -172,7 +172,7 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec JetTaggingUtils::get_ghostStatus(ghostFlavour ghostStruct){ +std::vector> JetTaggingUtils::get_jetconstituents(ghostFlavour ghostStruct){ + return ghostStruct.jet_constituents; +} + +ROOT::VecOps::RVec JetTaggingUtils::get_ghostStatus(ghostFlavour ghostStruct){ return ghostStruct.ghostStatus; } -std::vector JetTaggingUtils::get_MCindex(ghostFlavour ghostStruct){ +ROOT::VecOps::RVec JetTaggingUtils::get_MCindex(ghostFlavour ghostStruct){ return ghostStruct.MCindex; } diff --git a/analyzers/dataframe/JetTaggingUtils.h b/analyzers/dataframe/JetTaggingUtils.h index 20523bfbff..ae22fd377a 100644 --- a/analyzers/dataframe/JetTaggingUtils.h +++ b/analyzers/dataframe/JetTaggingUtils.h @@ -22,24 +22,33 @@ namespace JetTaggingUtils{ * Jet tagging interface utilities. */ - //Get flavour association of jet + /** Get flavour association of jet */ ROOT::VecOps::RVec get_flavour(ROOT::VecOps::RVec in, ROOT::VecOps::RVec MCin); - + + /** structure to output the following: + * - vector of parton and hadron flavours + * - FCCAnalyses jets clustered by user defined algorithm + * - jet constituent indices + * - ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) + * - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) + */ struct ghostFlavour { std::vector> flavour; JetClusteringUtils::FCCAnalysesJet jets; - std::vector ghostStatus; - std::vector MCindex; + std::vector> jet_constituents; + ROOT::VecOps::RVec ghostStatus; + ROOT::VecOps::RVec MCindex; }; - //Get ghost flavour (MC flavour) of jet described here: .. - //a struct is returned containing - // - vector of flavours: <, > - // - FCCAnalysesJet struct containing resulting jets (on which functions defined in JetClusteringUtils.cc can be called) - // - vector of jet constituent indices <, , ...> - // - ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) - // - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) + /** Get ghost flavour (MC flavour) of jet described here: .. + * a struct is returned containing + * - vector of flavours: <, > + * - FCCAnalysesJet struct containing resulting jets (on which functions defined in JetClusteringUtils.cc can be called) + * - vector of jet constituent indices <, , ...> + * - ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) + * - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) + */ struct get_ghostFlavour { get_ghostFlavour(int algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, float arg_add1=0, float arg_add2=0); @@ -54,13 +63,20 @@ namespace JetTaggingUtils{ ghostFlavour operator() (ROOT::VecOps::RVec Particle, ROOT::VecOps::RVec ind, std::vector pseudoJets, int partonFlag); }; + /** get the vector of flavours: <, > */ std::vector> get_flavour(ghostFlavour ghostStruct); + /** get the FCCAnalysesJet struct containing resulting jets (on which functions defined in JetClusteringUtils.cc can be called) */ JetClusteringUtils::FCCAnalysesJet get_jets(ghostFlavour ghostStruct); - std::vector get_ghostStatus(ghostFlavour ghostStruct); + /** get the vector of jet constituent indices <, , ...> */ + std::vector> get_jetconstituents(ghostFlavour ghostStruct); + + /** get ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) */ + ROOT::VecOps::RVec get_ghostStatus(ghostFlavour ghostStruct); - std::vector get_MCindex(ghostFlavour ghostStruct); + /** get the MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) */ + ROOT::VecOps::RVec get_MCindex(ghostFlavour ghostStruct); //Get b-tags with an efficiency applied ROOT::VecOps::RVec get_btag(ROOT::VecOps::RVec in, float efficiency, float mistag_c=0., float mistag_l=0., float mistag_g=0.); From baa08e45bb5a7e1ff67f581fdd33562d42976870 Mon Sep 17 00:00:00 2001 From: Kunal Gautam Date: Tue, 15 Mar 2022 17:06:25 +0100 Subject: [PATCH 5/6] struct updated --- analyzers/dataframe/JetTaggingUtils.cc | 12 ++++++---- analyzers/dataframe/JetTaggingUtils.h | 32 +++++++++++++++++--------- 2 files changed, 29 insertions(+), 15 deletions(-) diff --git a/analyzers/dataframe/JetTaggingUtils.cc b/analyzers/dataframe/JetTaggingUtils.cc index 031ec210b0..c8a4b685d3 100644 --- a/analyzers/dataframe/JetTaggingUtils.cc +++ b/analyzers/dataframe/JetTaggingUtils.cc @@ -62,10 +62,14 @@ JetTaggingUtils::get_flavour(ROOT::VecOps::RVec in, } -get_ghostFlavour::get_ghostFlavour(int algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, float arg_add1, float arg_add2) -{m_algo = algo; m_radius = arg_radius; m_exclusive = arg_exclusive; m_cut = arg_cut; m_sorted = arg_sorted; m_recombination = arg_recombination; m_add1 = arg_add1; m_add2 = arg_add2;} - -ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec Particle, ROOT::VecOps::RVec ind, std::vector pseudoJets, int partonFlag) { +get_ghostFlavour::get_ghostFlavour(int arg_algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, + float arg_add1, float arg_add2) +{m_algo = arg_algo; m_radius = arg_radius; m_exclusive = arg_exclusive; m_cut = arg_cut; m_sorted = arg_sorted; m_recombination = arg_recombination; m_add1 = arg_add1; m_add2 = arg_add2;} + +ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec Particle, + ROOT::VecOps::RVec ind, + std::vector pseudoJets, + int partonFlag) { diff --git a/analyzers/dataframe/JetTaggingUtils.h b/analyzers/dataframe/JetTaggingUtils.h index ae22fd377a..f131d44f27 100644 --- a/analyzers/dataframe/JetTaggingUtils.h +++ b/analyzers/dataframe/JetTaggingUtils.h @@ -50,17 +50,27 @@ namespace JetTaggingUtils{ * - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) */ struct get_ghostFlavour { - get_ghostFlavour(int algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, float arg_add1=0, float arg_add2=0); - - int m_algo = 0; ///< flag to select jet clustering algorithm defined in JetClustering.cc (0 = kt, 1 = antikt, 2 = cambridge, 3 = eekt, 4 = ee genkt, 5 = genkt, 6 = valencia, 7 = jade) - float m_radius = 0.5; ///< jet cone radius (note that this variable should be passed even when not used e.g. for eekt) - int m_exclusive = 0; ///< flag for exclusive jet clustering. Possible choices are 0=inclusive clustering, 1=exclusive clustering that would be obtained when running the algorithm with the given dcut, 2=exclusive clustering when the event is clustered (in the exclusive sense) to exactly njets, 3=exclusive clustering when the event is clustered (in the exclusive sense) up to exactly njets, 4=exclusive jets obtained at the given ycut - float m_cut = 5.; ///< pT cut for m_exclusive=0, dcut for m_exclusive=1, N jets for m_exlusive=2, N jets for m_exclusive=3, ycut for m_exclusive=4 - int m_sorted = 0; ///< pT ordering=0, E ordering=1 - int m_recombination = 0; ///< E_scheme=0, pt_scheme=1, pt2_scheme=2, Et_scheme=3, Et2_scheme=4, BIpt_scheme=5, BIpt2_scheme=6, E0_scheme=10, p_scheme=11 - float m_add1 = 0.; /// first additional parameter (to be used in selected clustering scheme) - float m_add2 = 0.; /// second additional parameter - ghostFlavour operator() (ROOT::VecOps::RVec Particle, ROOT::VecOps::RVec ind, std::vector pseudoJets, int partonFlag); + public: + get_ghostFlavour( int arg_algo = 0, float arg_radius = 0.5, int arg_exclusive = 0, float arg_cut = 5., int arg_sorted = 0, int arg_recombination = 0, + float arg_add1 = 0, float arg_add2 = 0 ); + + /// arg_algo: flag to select jet clustering algorithm defined in JetClustering.cc (0 = kt, 1 = antikt, 2 = cambridge, 3 = eekt, 4 = ee genkt, 5 = genkt, 6 = valencia, 7 = jade) + /// arg_radius: jet cone radius (note that this variable should be passed even when not used e.g. for eekt) + /// arg_exclusive: flag for exclusive jet clustering. Possible choices are 0=inclusive clustering, 1=exclusive clustering that would be obtained when running the algorithm with the given dcut, 2=exclusive clustering when the event is clustered (in the exclusive sense) to exactly njets, 3=exclusive clustering when the event is clustered (in the exclusive sense) up to exactly njets, 4=exclusive jets obtained at the given ycut + /// arg_cut: pT cut for m_exclusive=0, dcut for m_exclusive=1, N jets for m_exlusive=2, N jets for m_exclusive=3, ycut for m_exclusive=4 + /// arg_sorted: pT ordering=0, E ordering=1 + /// arg_recombination: E_scheme=0, pt_scheme=1, pt2_scheme=2, Et_scheme=3, Et2_scheme=4, BIpt_scheme=5, BIpt2_scheme=6, E0_scheme=10, p_scheme=11 + /// arg_add1: first additional parameter (to be used in selected clustering scheme) + /// arg_add2: second additional parameter + + ghostFlavour operator() (ROOT::VecOps::RVec Particle, + ROOT::VecOps::RVec ind, + std::vector pseudoJets, + int partonFlag); + + private: + int m_algo; float m_radius; int m_exclusive; float m_cut; int m_sorted; int m_recombination; + float m_add1; float m_add2; }; /** get the vector of flavours: <, > */ From f79255c571cd3eef0a5fcf5f9dea50b121f4f337 Mon Sep 17 00:00:00 2001 From: Eduardo Ploerer Date: Wed, 16 Mar 2022 18:08:48 +0100 Subject: [PATCH 6/6] Passing values by reference --- analyzers/dataframe/JetTaggingUtils.cc | 36 +++++++++++--------------- analyzers/dataframe/JetTaggingUtils.h | 34 ++++++++++++------------ 2 files changed, 33 insertions(+), 37 deletions(-) diff --git a/analyzers/dataframe/JetTaggingUtils.cc b/analyzers/dataframe/JetTaggingUtils.cc index c8a4b685d3..a0e00fdc12 100644 --- a/analyzers/dataframe/JetTaggingUtils.cc +++ b/analyzers/dataframe/JetTaggingUtils.cc @@ -62,22 +62,21 @@ JetTaggingUtils::get_flavour(ROOT::VecOps::RVec in, } -get_ghostFlavour::get_ghostFlavour(int arg_algo, float arg_radius, int arg_exclusive, float arg_cut, int arg_sorted, int arg_recombination, - float arg_add1, float arg_add2) +get_ghostFlavour::get_ghostFlavour(const int & arg_algo, const float & arg_radius, const int & arg_exclusive, const float & arg_cut, const int & arg_sorted, const int & arg_recombination, + const float & arg_add1, const float & arg_add2) {m_algo = arg_algo; m_radius = arg_radius; m_exclusive = arg_exclusive; m_cut = arg_cut; m_sorted = arg_sorted; m_recombination = arg_recombination; m_add1 = arg_add1; m_add2 = arg_add2;} -ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec Particle, - ROOT::VecOps::RVec ind, - std::vector pseudoJets, - int partonFlag) { - +ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (const ROOT::VecOps::RVec & Particle, + const ROOT::VecOps::RVec & ind, + std::vector & pseudoJets, + const int & partonFlag) { ghostFlavour result; unsigned int index = pseudoJets.size(); - ROOT::VecOps::RVec pdg(pseudoJets.size(),0); - ROOT::VecOps::RVec ghostStatus(pseudoJets.size(),0); + ROOT::VecOps::RVec pdg(pseudoJets.size(),0); + ROOT::VecOps::RVec ghostStatus(pseudoJets.size(),0); ROOT::VecOps::RVec MCindex(pseudoJets.size(),-1); @@ -176,24 +175,23 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec> flav_vector; + std::vector> flav_vector; - std::vector partonFlavs; - std::vector hadronFlavs; + std::vector partonFlavs; + std::vector hadronFlavs; for (auto& consti_index : jetconstituents) { - float partonFlav = 0; + int partonFlav = 0; float partonMom2 = 0; float partonMom2_b = 0; float partonMom2_c = 0; - float hadronFlav = 0; + int hadronFlav = 0; float hadronMom2_b = 0; float hadronMom2_c = 0; @@ -241,7 +239,6 @@ ghostFlavour JetTaggingUtils::get_ghostFlavour::operator() (ROOT::VecOps::RVec> JetTaggingUtils::get_flavour(ghostFlavour ghostStruct){ +std::vector> JetTaggingUtils::get_flavour(ghostFlavour ghostStruct){ return ghostStruct.flavour; } @@ -258,11 +255,8 @@ JetClusteringUtils::FCCAnalysesJet JetTaggingUtils::get_jets(ghostFlavour ghostS return ghostStruct.jets; } -std::vector> JetTaggingUtils::get_jetconstituents(ghostFlavour ghostStruct){ - return ghostStruct.jet_constituents; -} -ROOT::VecOps::RVec JetTaggingUtils::get_ghostStatus(ghostFlavour ghostStruct){ +ROOT::VecOps::RVec JetTaggingUtils::get_ghostStatus(ghostFlavour ghostStruct){ return ghostStruct.ghostStatus; } diff --git a/analyzers/dataframe/JetTaggingUtils.h b/analyzers/dataframe/JetTaggingUtils.h index f131d44f27..24a3ed7d74 100644 --- a/analyzers/dataframe/JetTaggingUtils.h +++ b/analyzers/dataframe/JetTaggingUtils.h @@ -26,17 +26,16 @@ namespace JetTaggingUtils{ ROOT::VecOps::RVec get_flavour(ROOT::VecOps::RVec in, ROOT::VecOps::RVec MCin); /** structure to output the following: - * - vector of parton and hadron flavours + * - vector of parton and hadron flavours (Note that this is not an RVec but rather std vec, since there were issues with RDataFrame Snapshot) * - FCCAnalyses jets clustered by user defined algorithm * - jet constituent indices * - ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) * - MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) */ struct ghostFlavour { - std::vector> flavour; + std::vector> flavour; JetClusteringUtils::FCCAnalysesJet jets; - std::vector> jet_constituents; - ROOT::VecOps::RVec ghostStatus; + ROOT::VecOps::RVec ghostStatus; ROOT::VecOps::RVec MCindex; }; @@ -51,8 +50,14 @@ namespace JetTaggingUtils{ */ struct get_ghostFlavour { public: - get_ghostFlavour( int arg_algo = 0, float arg_radius = 0.5, int arg_exclusive = 0, float arg_cut = 5., int arg_sorted = 0, int arg_recombination = 0, - float arg_add1 = 0, float arg_add2 = 0 ); + get_ghostFlavour(const int & arg_algo = 0, + const float & arg_radius = 0.5, + const int & arg_exclusive = 0, + const float & arg_cut = 5., + const int & arg_sorted = 0, + const int & arg_recombination = 0, + const float & arg_add1 = 0, + const float & arg_add2 = 0 ); /// arg_algo: flag to select jet clustering algorithm defined in JetClustering.cc (0 = kt, 1 = antikt, 2 = cambridge, 3 = eekt, 4 = ee genkt, 5 = genkt, 6 = valencia, 7 = jade) /// arg_radius: jet cone radius (note that this variable should be passed even when not used e.g. for eekt) @@ -62,11 +67,11 @@ namespace JetTaggingUtils{ /// arg_recombination: E_scheme=0, pt_scheme=1, pt2_scheme=2, Et_scheme=3, Et2_scheme=4, BIpt_scheme=5, BIpt2_scheme=6, E0_scheme=10, p_scheme=11 /// arg_add1: first additional parameter (to be used in selected clustering scheme) /// arg_add2: second additional parameter - - ghostFlavour operator() (ROOT::VecOps::RVec Particle, - ROOT::VecOps::RVec ind, - std::vector pseudoJets, - int partonFlag); + + ghostFlavour operator() (const ROOT::VecOps::RVec & Particle, + const ROOT::VecOps::RVec & ind, + std::vector & pseudoJets, + const int & partonFlag); private: int m_algo; float m_radius; int m_exclusive; float m_cut; int m_sorted; int m_recombination; @@ -74,16 +79,13 @@ namespace JetTaggingUtils{ }; /** get the vector of flavours: <, > */ - std::vector> get_flavour(ghostFlavour ghostStruct); + std::vector> get_flavour(ghostFlavour ghostStruct); /** get the FCCAnalysesJet struct containing resulting jets (on which functions defined in JetClusteringUtils.cc can be called) */ JetClusteringUtils::FCCAnalysesJet get_jets(ghostFlavour ghostStruct); - - /** get the vector of jet constituent indices <, , ...> */ - std::vector> get_jetconstituents(ghostFlavour ghostStruct); /** get ghostStatus (0 = not ghost, 1 = ghost parton, 2 = ghost hadron) */ - ROOT::VecOps::RVec get_ghostStatus(ghostFlavour ghostStruct); + ROOT::VecOps::RVec get_ghostStatus(ghostFlavour ghostStruct); /** get the MCindex vector of indices of ghosts to MC collection (and -1 if not a ghost) */ ROOT::VecOps::RVec get_MCindex(ghostFlavour ghostStruct);