Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHPy6JetTrigger.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHPy6JetTrigger.cc
1 #include "PHPy6JetTrigger.h"
2 #include "PHPy6GenTrigger.h"
3 
4 #pragma GCC diagnostic push
5 #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
6 #include <HepMC/GenEvent.h>
7 #pragma GCC diagnostic pop
8 
9 #include <HepMC/GenParticle.h> // for GenParticle
10 #include <HepMC/SimpleVector.h> // for FourVector
11 
12 // fastjet includes
13 #include <fastjet/ClusterSequence.hh>
14 #include <fastjet/JetDefinition.hh>
15 #include <fastjet/PseudoJet.hh>
16 
17 #include <cmath> // for sqrt
18 #include <cstdlib> // for abs
19 #include <iostream> // for operator<<, endl, basic_ostream
20 #include <memory> // for allocator_traits<>::value_type
21 #include <utility> // for swap
22 #include <vector> // for vector
23 
24 using namespace std;
25 
26 //__________________________________________________________
28  : PHPy6GenTrigger(name)
29 {
30 }
31 
33 {
34  if (Verbosity() > 0) PrintConfig();
35 }
36 
37 bool PHPy6JetTrigger::Apply(const HepMC::GenEvent *evt)
38 {
39  // Loop over all particles in the event
40  int idx = 0;
41  std::vector<fastjet::PseudoJet> pseudojets;
42  for (HepMC::GenEvent::particle_const_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p)
43  {
44  idx++;
45  if (((*p)->status() != 1) != 0) continue;
46 
47  // remove some particles (muons, taus, neutrinos)...
48  // 12 == nu_e
49  // 13 == muons
50  // 14 == nu_mu
51  // 15 == taus
52  // 16 == nu_tau
53  if ((abs(((*p)->pdg_id())) >= 12) && (abs(((*p)->pdg_id())) <= 16)) continue;
54 
55  // acceptance... _etamin,_etamax
56  if (((*p)->momentum().px() == 0.0) && ((*p)->momentum().py() == 0.0)) continue; // avoid pt=0
57  if ((((*p)->momentum().pseudoRapidity()) < m_theEtaLow) ||
58  (((*p)->momentum().pseudoRapidity()) > m_theEtaHigh)) continue;
59 
60  fastjet::PseudoJet pseudojet((*p)->momentum().px(),
61  (*p)->momentum().py(),
62  (*p)->momentum().pz(),
63  (*p)->momentum().e());
64  pseudojet.set_user_index(idx);
65  pseudojets.push_back(pseudojet);
66  }
67 
68  // Call FastJet
69 
70  fastjet::JetDefinition *jetdef = new fastjet::JetDefinition(fastjet::antikt_algorithm, m_R, fastjet::E_scheme, fastjet::Best);
71  fastjet::ClusterSequence jetFinder(pseudojets, *jetdef);
72  std::vector<fastjet::PseudoJet> fastjets = jetFinder.inclusive_jets();
73  delete jetdef;
74 
75  bool jetFound = false;
76  double max_pt = -1;
77  for (unsigned int ijet = 0; ijet < fastjets.size(); ++ijet)
78  {
79  const double pt = sqrt(pow(fastjets[ijet].px(), 2) + pow(fastjets[ijet].py(), 2));
80 
81  if (pt > max_pt) max_pt = pt;
82 
83  vector<fastjet::PseudoJet> constituents = fastjets[ijet].constituents();
84  const int nconst = constituents.size();
85 
86  if (pt > m_minPt && nconst >= m_nconst)
87  {
88  jetFound = true;
89  break;
90  }
91  }
92 
93  if (Verbosity() > 2)
94  {
95  cout << "PHPy6JetTrigger::Apply - max_pt = " << max_pt << ", and jetFound = " << jetFound << endl;
96  }
97 
98  return jetFound;
99 }
100 
101 void PHPy6JetTrigger::SetEtaHighLow(double etaHigh, double etaLow)
102 {
103  m_theEtaHigh = etaHigh;
104  m_theEtaLow = etaLow;
105 
107  {
109  }
110 }
111 
113 {
114  cout << "---------------- PHPy6JetTrigger::PrintConfig --------------------" << endl;
115 
116  cout << " Particles EtaCut: " << m_theEtaLow << " < eta < " << m_theEtaHigh << endl;
117  cout << " Minimum Jet pT: " << m_minPt << " GeV/c" << endl;
118  cout << " Anti-kT Radius: " << m_R << endl;
119  cout << "-----------------------------------------------------------------------" << endl;
120 }