2 #include "include/smearing.h"
8 void SortEvent_ATLAS_GAMBIT(Pythia8::Event &event,
HEP::Event &OutEvt,
int npileups, std::vector<HEP::PileupEvent*> &pileups,
double &Rjet, std::mt19937 &engine){
10 fastjet::JetDefinition jetDef(fastjet::antikt_algorithm, Rjet);
11 std::vector <HEP::Jet*> outjets;
13 std::vector <fastjet::PseudoJet> inclusiveJets, jets;
15 std::vector <fastjet::PseudoJet> pileupJets, pJInputs;
20 std::vector <fastjet::PseudoJet> fjInputs, bquark_momenta;
23 bquark_momenta.resize(0);
28 for(
int mpart = event.size()-1; mpart>=0 ; mpart--)
31 Pythia8::Particle* prt = &
event[mpart];
43 if(prt->isHadron() && prt->isVisible())
45 Pythia8::Vec4 pos = prt->vProd();
47 if (pos.pAbs() > 0.0001)
49 newpart->set_prod(
HEP::P4(pos.px(),pos.py(),pos.pz(),pos.e()));
56 newpart->set_3Q(prt->chargeType());
57 OutEvt.add_charged_hadron(newpart);
69 OutEvt.add_neutral_hadron(newpart->
mom());
81 Pythia8::Vec4 pos = prt->vProd();
84 if (pos.pAbs() > 0.0001)
86 newpart->set_prod(
HEP::P4(pos.px(),pos.py(),pos.pz(),pos.e()));
92 if(apid == 11 || apid == 13 || apid == 22 )
94 if(!fromHadron(mpart,event))
98 if(apid == 11 || apid ==13 )
100 if(pid > 0) {newpart->set_3Q(-3);}
else {newpart->set_3Q(3);}
113 OutEvt.add_neutral_hadron(newpart->
mom());
118 OutEvt.add_charged_hadron(newpart);
131 newpart->set_3Q(prt->chargeType());
143 if( pid == 12 || pid == 14 || pid == 16 || pid == 18 )
147 if((apid > 1000000) && (apid < 9000000))
152 Pythia8::Vec4 v4=prt->p();
153 fastjet::PseudoJet particleTemp(v4.px(),v4.py(),v4.pz(),v4.e());
154 fjInputs.push_back(particleTemp);
188 std::vector<int> daughters = prt->daughterList();
189 for(
int di = 0; di < daughters.size(); di ++)
191 if (event[daughters[di]].idAbs() == 15)
203 Pythia8::Vec4 pos = prt->vProd();
204 if (pos.pAbs() > 0.001)
206 newpart->set_prod(
HEP::P4(pos.px(),pos.py(),pos.pz(),pos.e()));
209 if (pos.pAbs() > 0.001)
211 newpart->set_decay(
HEP::P4(pos.px(),pos.py(),pos.pz(),pos.e()));
257 std::vector<int> daughters = prt->daughterList();
258 for(
int di = 0; di < daughters.size(); di ++)
260 if (event[daughters[di]].idAbs() == 5)
270 Pythia8::Vec4 v4 = prt->p();
272 fastjet::PseudoJet tpjet(v4.px(),v4.py(),v4.pz(),v4.e());
274 bquark_momenta.push_back(tpjet);
284 std::uniform_int_distribution<> rd_el(0,pileups.size()-1);
286 static std::normal_distribution<double> zdist(0,53.0);
287 static std::normal_distribution<double> tdist(0,160e-12);
296 fastjet::ClusterSequence clustSeq(fjInputs, jetDef);
297 inclusiveJets = clustSeq.inclusive_jets(1.0);
298 jets = sorted_by_pt(inclusiveJets);
313 for(
int jp=0; jp< npileups; jp++)
315 int pileup_number=rd_el(engine);
317 double dz=zdist(engine);
318 double dt=tdist(engine);
321 OutEvt.
add_particles(pileups[pileup_number]->translate_particles(dz,dt));
322 OutEvt.add_charged_hadrons(pileups[pileup_number]->translate_charged_hadrons(dz,dt));
325 OutEvt.add_neutral_hadrons(pileups[pileup_number]->get_neutral_hadrons());
395 for(
auto const& jjet: jets)
398 for(
auto const& bquark : bquark_momenta)
400 if(bquark.delta_R(jjet) < Rjet) { btag =
true;
break;}
409 outjets.push_back(tjet);
433 void filleventPYTHIA_ATLAS_GAMBIT(
HEP::Event &OutEvt, Pythia8::Event &evt,
double weight,
SettingsContainer &runsettings,std::vector<HEP::PileupEvent*> &pileups, std::mt19937 &engine)
439 OutEvt.set_weight(weight);
440 int npileups = runsettings.npileup;
441 double dR = runsettings.Rjet;
445 SortEvent_ATLAS_GAMBIT(evt,OutEvt,npileups,pileups,dR,engine);
465 simplecalcMissingMomentum(OutEvt,engine);
Simple event class, separating particles into classes.
Definition: heputil.h:227
void add_particles(const std::vector< Particle * > &ps)
Definition: heputil.h:584
const std::vector< Particle * > & muons() const
Get prompt muons.
Definition: heputil.h:661
const std::vector< Particle * > & electrons() const
Get prompt electrons.
Definition: heputil.h:651
void add_particle(Particle *p)
Definition: heputil.h:523
void set_jets(const std::vector< Jet * > &jets)
Set the jets collection.
Definition: heputil.h:709
A 4-momentum class for vectors.
Definition: Vectors.h:45
Definition: Particle.h:24
const P4 & mom() const
Get the 4 vector.
Definition: Particle.h:109
void set_prompt(bool isprompt=true)
Set promptness.
Definition: Particle.h:178
Class to hold settings, read from the YAML file, such as number of cores, deltaR for the jets,...
Definition: running.h:65
void smearJets(std::vector< HEP::Jet * > &jets, std::mt19937 &engine)
Randomly smear the supplied jets' momenta by parameterised resolutions.
Definition: smearing.h:491
void smearElectronEnergy(std::vector< HEP::Particle * > &electrons, std::mt19937 &engine)
Randomly smear the supplied electrons' momenta by parameterised resolutions.
Definition: smearing.h:413
void smearMuonMomentum(std::vector< HEP::Particle * > &muons, std::mt19937 &engine)
Randomly smear the supplied muons' momenta by parameterised resolutions.
Definition: smearing.h:460